Free Access
Issue
A&A
Volume 642, October 2020
Article Number A129
Number of page(s) 9
Section Extragalactic astronomy
DOI https://doi.org/10.1051/0004-6361/202038052
Published online 13 October 2020

© ESO 2020

1. Introduction

Rapid variability has been one of the identifying criteria of active galactic nuclei (AGNs), which are powered by accretion on the central super massive black hole (SMBH). In particular the subset of AGNs that are designated as radio-loud have been found to exhibit a diverse range of observational behaviors across all the electromagnetic bands. The blazar subclass of radio-loud AGNs host powerful large-scale relativistic jets of plasma that is nearly pointed toward us (Urry & Padovani 1995). They emit an entirely jet-dominated radiation, spread across the entire accessible electromagnetic (EM) spectrum (e.g., Abdo et al. 2010). Blazars are known for their enormous dynamic range of rarely repeating observational behaviors. In the extragalactic sky, they are the most prominent and persistent broadband nonthermal emitters. Observational studies have revealed them to be variable on all time scales, all the way from the shortest time scale allowed by the observing facilities to the longest one allowed by the available data records (e.g., Dey et al. 2018; Shukla et al. 2018, and references therein).

In general, accretion-powered sources exhibit variability on a broad range of time scales and exhibit a diverse range of variability behaviors. Despite this, their variability has been found to share some common statistical variability properties (Scaringi et al. 2015). Given the fact that accretion-powered sources encompass both compact and noncompact astrophysical objects of all mass scales from proto-stars to AGNs, with intrinsically very different physical conditions, processes, and mechanisms, the inferred similarity has been argued to be related to the accretion physics. Blazars’ broadband flux variability, despite being dominated by jet emission, also broadly follow these properties (Kushwaha et al. 2016, 2017; Bhatta & Dhital 2020).

Blazar multi-wavelength (MW) flux variability is, in general, stochastic, but quasi-periodic variations (QPOs) have been reported occasionally (Valtaoja et al. 1985; Carrasco et al. 1985; Quirrenbach et al. 1991; Urry et al. 1993; Sillanpaa et al. 1996; Heidt & Wagner 1996; Raiteri et al. 2003; Espaillat et al. 2008; Gupta et al. 2008, 2009, 2019; Lachowicz et al. 2009; King et al. 2013; Sandrinelli et al. 2014, 2016, 2017; Gupta 2014, 2018; Ackermann et al. 2015; Bhatta 2017, 2019; Zhang et al. 2017; Hong et al. 2018, and references therein). The suggested QPO time scales also have a huge range, similar to the variability time scales shown by blazars, that is, from a few tens of minutes to hours, days, months and even years. (e.g., Kellermann & Pauliny-Toth 1968; Stein et al. 1976; Miller et al. 1989; Edelson et al. 1991, 1995; Sembay et al. 1993; Brinkmann et al. 1994; Wagner & Witzel 1995; Sillanpaa et al. 1996; Heidt & Wagner 1996; Gupta et al. 2008; Gaur et al. 2012, 2015; Dey et al. 2018; Kushwaha 2020, and references therein). In addition to QPOs, as observed in Galactic black hole X-ray binaries, that are normally attributed to accretion, some QPOs in blazars are also expected from the cosmological hierarchical structure formation, suggesting the formation of binary supermassive black holes. However, this latter type is expected to be persistent, compared to the accretion induced QPOs which are generally transient. Further, given the complex physics and not yet fully understood connections between different constituents and jets of AGNs, their imprints, if any, are expected to be reflected in the jet emission. In this context, the continuous monitoring capability of the Fermi-LAT instrument provides an excellent facility to explore bright extragalactic γ-ray sources, most of which are blazars, on both short and long time scales (e.g., Shukla et al. 2018; Kushwaha et al. 2017; Gupta et al. 2019).

CTA 102 is a flat spectrum radio quasar (FSRQ) – the blazar subclass characterized by prominent broad emission lines – located at the redshift of 1.032 (Schmidt 1965). It was first discovered in the radio survey by Owens Valley Radio Observatory (OVRO) at 960 MHz (serial number 102; Harris & Roberts 1960) and its optical counterpart was identified by Sandage & Wyndham (1965). Like other blazars and radio-loud AGNs, it has been explored across multiple EM bands and found to exhibit the characteristics of radio-loud AGNs, for example, compact star-like appearances in images (Sandage & Wyndham 1965), rapid and high broadband variability (Pica et al. 1988; Meyer et al. 2019; D’Ammando et al. 2019), a high (> 3%) and variable optical polarization degree (Moore & Stockman 1981; Casadio et al. 2015), high brightness temperature (Fromm et al. 2013) and nonthermal spectrum (and references therein Gasparyan et al. 2018; Raiteri et al. 2017) except that it is more luminous and has broad emission lines in the optical band only (Larionov et al. 2016). Optical observations before the year 2000 suggest moderate variability on long terms (Pica et al. 1988), but stronger optical variability has been reported after this (Osterman Meyer et al. 2009). From mid-2011, the source entered an extended high activity period, exhibiting its highest ever reported flux state from the end of 2016 to the beginning of 2017 (e.g., Raiteri et al. 2017; D’Ammando et al. 2019)

At γ-ray energies, CTA 102 has been detected by the EGRET (Hartman et al. 1999) and COMPTEL (Blom et al. 1995) telescopes onboard the Compton Gamma Ray Observatory (CGRO). It also featured in the first three months of Fermi bright sources list (4C +11.69; Abdo et al. 2010) CTA 102, however, became very active in the Fermi-LAT band only in 2012, reaching a flux level detectable on daily time scales. Since then its γ-ray as well as MW activity kept increasing, reaching a peak flux of ∼2 × 10−5 ph cm−2 s−1 in the LAT band at the end of 2017 (Shukla et al. 2018; Zacharias et al. 2019; Meyer et al. 2019; D’Ammando et al. 2019). It was then the brightest object in the γ-ray sky with an isotropic luminosity of ∼3.25 × 1050 erg s−1 (Gasparyan et al. 2018) on daily time scales. The very high flux allowed a detailed study of variability down to time scales as small as a few minutes (Shukla et al. 2018). A MW spectral study by Zacharias et al. (2019) argue, that the high activity spread over a few months is due to ablation of a cloud by the relativistic jet and the resulting hadronic interactions. Another possible explanation is that the high MW activity is likely to involve the interaction between a superluminal component and stationary core ∼0.1 mas from the core, considered to be a recollimation shock (Casadio et al. 2019). On the other hand, a detailed and systematic study by Raiteri et al. (2017) shows that the increased MW activity provides hardly any evidence of spectral evolution and thus they explained its long term spectral and temporal variability as primarily resulting from variations in the Doppler factor of an inhomogeneous jet.

In this work, we report the detection of a QPO, in both the γ−ray and optical R-band light curves of the FSRQ CTA 102, during its highest reported activity period from the end of 2016 to the beginning of 2017 (MJD: 57710–57790). We found a highly significant QPO of ∼7.6 days. We look at the data reduction methodology in Sect. 2. The analysis procedure and the results are in Sects. 3 and 4. We discuss the possible physical processes explaining the QPO in Sect. 5 and present our conclusions in Sect. 6.

2. Data acquisition

2.1. Gamma-ray data

The γ-ray data used in our work, is taken from the Large Area Telescope (LAT) facility, onboard the Fermi observatory. LAT is an imaging telescope working on the principle of pair conversion and is sensitive to photon energies > 20 MeV. It has a large angular field of view of 2 sr and covers the entire sky every ∼90 min (Atwood et al. 2009).

We used the PASS8 (P8R2) processed events data of CTA 102 between MJD 57710–57790. The PASS8 data is an improved reconstruction of the entire LAT events and provides a significant improvement in the data quality (Fermi-LAT Collaboration 2020). As per the recommendation of the instrument team, we considered only “SOURCE” class events (evclass = 128, evtype = 3) within the energy range 0.1−300 GeV from a circular region of 15° centered on the source. A zenith angle cutoff of 90° was applied during the event selection to avoid contamination from the Earth’s limb. Following this, the good time intervals were generated using the standard expression “(DATA_QUAL>0)&&(LAT_CONFIG==1)”. Subsequently, an exposure map was calculated with angular coverage 10° larger than the region of interest. An XML file containing the spectral shape of point sources, in the chosen field of view, from the third Fermi-LAT source catalog (3FGL; Acero et al. 2015), was generated. Additionally, it also included the Galactic and extra-galactic contribution through their respective template files “gll_iem_v06” and “iso_P8R2_SOURCE_V6_v06”. Finally, the optimization over the input XML spectral file was performed using the unbinned likelihood method which comes bundled with the LAT analysis software.

The γ-ray light curve of the source was extracted following the above-mentioned procedure for each time interval. First, we extracted the daily light curve of the source by iteratively removing nonsignificant sources during the likelihood fit as described in Kushwaha et al. (2014) until the analysis converged. We then used the daily best-fit source XML files to extract the 3h light curve. Following the data reduction techniques mentioned above, the γ-ray analysis was performed in sequential time steps, with time bins of 3 hours to obtain the light curve of the source during its flaring state. Bins of three-hour length, optimized both cadence and coverage for the purpose of our analysis; they were long enough to give adequate signal-to-noise ratios (test statistics > 9) in nearly all bins, and short enough to give a large enough number of points for good light curve modeling. When rebinned to 24 hour segments, our light curve was in very good agreement with the Fermi-LAT collaboration daily quick-look product. As discussed below, we checked that the key results were essentially identical when somewhat longer bins were considered.

2.2. Optical (R-band) data

The R passband optical photometric data presented here is collected from the various telescopes around the globe, under the whole earth blazar telescope (WEBT) observational monitoring campaigns of the blazar CTA 102. The original R passband optical photometric data of these observational campaigns have been published (Raiteri et al. 2017; D’Ammando et al. 2019), and were kindly provided to us by C. M. Raiteri. The details concerning the R-band optical photometric data and analysis methods have been published (Raiteri et al. 2017; D’Ammando et al. 2019). We received the data in magnitude vs time format which we converted to flux vs time. For consistency, the optical fluxes were also binned into 3-hour segments.

3. Analysis

The initial hint of quasi periodicity came from the visual inspection of the optical and γ-ray light curve of CTA 102, during its 2016-17 flare (Fig. 1). We employed several independent techniques to identify and quantify periodicity. The individual techniques are elaborated below.

thumbnail Fig. 1.

CTA 102 light curves during the 2016–17 flare. (a) Fermi-LAT (0.1−300 GeV) light curve. The black line gives the optimal block representation of the light curve. A peak is considered to be significant if the block height is greater than a set threshold. (b) Optical R band light curve and its optimal block representation. The QPO peaks are denoted by red arrows. (c) γ-ray and R-band emission DCF. (d) Optical polarization angle during the detected QPO.

3.1. Light curve inspection

Visual inspection of the light curves hinted at quasi-periodicity, owing to several apparently equispaced features (peaks) in the light curve. However, a naïve visual inspection is often subject to bias, hence more rigorous methods are necessary to analyze features in the light curve. To properly define the features, we use the Bayesian Block (BB) representation of the light curve (Scargle et al. 2013). This involves segregating the light curve into its optimal block representation, by maximizing the fitness function among all segmentations of the light curve. This assumes the photon number distribution to be Poissonian and the False Alarm Probability (FAP) of the representation was <0.01. The block representation of the light curve helps in identifying the location and extent of significant peaks in the light curve. We define a peak as, a block (B) where the block height (average flux) is greater than both the preceding and succeeding blocks. We also obtain the extent of the peak by traveling downward from the peak in both forward and backward directions, and the peak flux (ℱp) from the height of the block B. By this definition, peaks can occur even in a relatively low flux state. Thus, we only statistically associate a peak with a flare, if ℱp is greater than 2σ from the mean of the global light curve. Meyer et al. (2019) gives the global properties of the γ-ray light curve and further details on the peak finding method. For the present analysis, we consider a γ-ray peak to be statistically significant, if ℱp > 6.6 × 10−6 ph cm−2 s−1. After identifying the significant peaks in the light curves, the relative separation between the peaks can help detect quasi-periodicity.

This analysis methodology for peaks does not trivially translate to the R-band data, since the long-term properties of this light curve (mean and variance) are unknown, even though some estimates are available at earlier times, for the B and V bands (Véron-Cetty & Véron 2010) as well as the G, H and K bands (Cutri et al. 2003). To alleviate this issue, we used a Discrete Correlation Function (DCF; Edelson & Krolik 1988) to quantify the correlation between the γ-ray and R-band optical emission. The presence of a significant correlation at zero lag would indicate quasi-periodicity in the R-band emission if the γ-ray emission shows a QPO (or vice versa). To calculate the DCF, we first subtract a linear baseline (Welsh 1999) and then compute the Unbinned DCF (UDCF) by

UDCF ij = ( F γ ( t i ) F γ ¯ ) ( F R ( t j ) F R ¯ ) σ 2 ( F γ ) σ 2 ( F R ) , $$ \begin{aligned} \mathrm{UDCF}_{i j}=\frac{\left(\mathcal{F} _{\gamma }\left(t_{i}\right)-\overline{\mathcal{F} _{\gamma }}\right)\left(\mathcal{F} _{R}\left(t_{j}\right)-\overline{\mathcal{F} _{R}}\right)}{\sqrt{\sigma ^{2}\left(\mathcal{F} _{\gamma }\right) \sigma ^{2}\left(\mathcal{F} _{R}\right)}}, \end{aligned} $$(1)

here ℱA(ti) is the A band flux at time ti and F ¯ $ \overline{\mathcal{F} } $ and σ2(ℱ) are its mean and variance, respectively. We obtained the DCF by binning the points in the range τ − Δτ/2 ≤ Δtij ≤ τ + Δτ/2. Here τ is the time lag, Δτ is the bin width and Δtij = titj, giving

D C F ( τ ) = 1 n k = 1 n U D C F k ( τ ) ± 1 n 1 k = 1 n ( U D C F k D C F ( τ ) ) 2 , $$ \begin{aligned} \begin{aligned}&\mathrm{D C F}(\tau ) =\frac{1}{n} \sum _{k=1}^n \mathrm{U D C F}_{k}(\tau ) \\&\qquad \qquad \pm \frac{1}{n-1} \sqrt{\sum _{k=1}^{n}\left(\mathrm{U D C F}_{k}-\mathrm{D C F}(\tau )\right)^{2}}, \end{aligned} \end{aligned} $$(2)

where n is the number of points binned following the binning criteria. Since blazar emissions are variable, we calculate the means and variances in Eq. (1) using only the points falling within given lag bins (White & Peterson 1994).

3.2. Power spectral analysis

Light curve inspection, including BB modeling, can only provide an idea about periodicity in a light curve. Proper quantification of a possible QPO is necessary for making any claim regarding the emission mechanisms at play. Power spectral density (PSD) or periodograms were used to compute the power emitted at different frequencies. PSDs are mod-squares of the Discrete Fourier transform (DFT) of the source light curve. If emission power is significantly higher in a particular frequency, then we have evidence that there may be a periodic component to the emission. The DFT algorithm needs uniformly sampled data points (flux in the present case), which is very unlikely for any astrophysical source. We instead use a Lomb-Scargle Periodogram (LSP, Lomb 1976; Scargle 1982), to take into account irregular sampling. This involves, fitting the light curve with sinusoidal functions, with different frequencies, and constructing a periodogram from the goodness of the fit. This method can directly detect persistent periodicities; however (as in this work) astrophysical periodicities could be transient. Transient periodicities could arise from effects in the blazar jet or accretion disk that are short-lived. A LSP is insensitive to transient periodicities unless we already know the extent of the periodicity stretch a priori.

To detect and quantify any transient QPOs, we used the Weighted Wavelet Z-transform method (WWZ, Foster 1996). This method decomposes the data into time and frequency domains (WWZ maps), by convolving the light curve with a time and frequency-dependent kernel. The present analysis uses the Morlet Kernel (Grossmann & Morlet 1984) with analytical form f(ω[tτ]) = exp((tτ) − 2(tτ)2). The WWZ map is then given by:

W ( ω , τ ; x ( t ) ) = ω 1 / 2 x ( t ) f ( ω ( t τ ) ) d t = ω 1 / 2 x ( ω 1 z + τ ) f ( z ) d z . $$ \begin{aligned} \begin{aligned}&W(\omega , \tau ; x(t)) = \omega ^{1 / 2} \int x(t) f^{*}(\omega (t-\tau )) \mathrm{d} t \\&\qquad \qquad \quad =\omega ^{-1 / 2} \int x\left(\omega ^{-1} z+\tau \right) f^{*}(z) \mathrm{d} z . \end{aligned} \end{aligned} $$(3)

Here f is the complex conjugate of the wavelet kernel f, ω is the scale factor (frequency), τ is the time-shift, and z = ω(tτ). The kernel acts as a windowed DFT, with the window exp(−2(tτ)2), with the window size depending on ω and the constant c. In the present analysis, we obtain the best value of the parameter c by matching the time-averaged WWZ and the LSP in the same time period. A WWZ map has the advantage of locating both, any dominant periods and their spans in time and hence is useful in locating the transient periodicities.

3.3. Light curve modeling

Blazar light curves (like any time-series) can be modeled using stochastic models. Modeling the light curve can reveal periodicity in the emission process. When the present emission depends on the past emissions (yielding fewer drastic jumps in the emission), the model used is an Autoregressive (AR) model. It is mathematically expressed as F ( t i ) = j = 1 p θ j F ( t i j ) + ϵ ( t i ) $ \mathcal{F} (t_{i}) = \sum\nolimits _{j=1}^p \theta _j\ \mathcal{F} (t_{i-j}) + \epsilon (t_{i}) $, where ℱ(ti) is the emission at time ti which depends on the p prior emissions, with θi the auto-regression coefficients, and ϵ(t) is a normally distributed random variable (fluctuation or forecasting error). Similarly, when the present emission depends on the past fluctuations, the model used is a Moving Average (MA) model, where F ( t i ) = j = 1 q ϕ j ϵ ( t i j ) + ϵ ( t i ) $ \mathcal{F} (t_{i}) = \sum\nolimits _{j=1}^q \phi _j \epsilon (t_{i-j}) + \epsilon (t_{i}) $ and as before, ϵ(ti) are normally distributed random variables.

The combined Autoregressive Moving Average (ARMA) model can nicely explain stationary timeseries; however, blazar light curves are often nonstationary. To deal with this issue, it is helpful to transform the light curve by successive differencing (Δd), defined as: Δℱ(ti) = ℱ(ti) − ℱ(ti−1), prior to modeling. This more complex model is then termed an Autoregressive (AR) Integrated (I) Moving Average (MA) or ARIMA(p, d, q) model (Scargle 1981; Feigelson et al. 2018), where the “integrated” part in the model accounts for the successive differencing. It is given by

Δ d F ( t i ) = j = 1 p θ j Δ d F ( t i j ) + j = 1 q ϕ j ϵ ( t i j ) + ϵ ( t i ) or , ( 1 j = 1 p θ j L j ) Δ d F ( t i ) = ( 1 j = 1 q ϕ j L j ) ϵ ( t i ) , $$ \begin{aligned} \begin{aligned}&\Delta ^d \mathcal{F} (t_{i}) = \sum _{j=1}^p \theta _j\ \Delta ^d\mathcal{F} (t_{i-j}) + \sum _{j=1}^q \phi _j \epsilon (t_{i-j}) + \epsilon (t_i)\\&\mathrm{or},\ \left(1 - \sum _{j=1}^p \theta _j \mathcal{L} ^j \right)\Delta ^d \mathcal{F} (t_i) = \left(1 - \sum _{j=1}^q \phi _j \mathcal{L} ^j \right) \epsilon (t_i) , \end{aligned} \end{aligned} $$(4)

where the successive differencing operator is defined as: Δd = (1 − ℐ)d. In the above equation, the second representation of the ARIMA model uses a lag operator defined as ℒkℱ(ti) = ℱ(tik) and ℒk ϵ(ti) = ϵ(tik) and p, q and d are the order of AR, MA, and differencing respectively.

An ARIMA(p,d,q) approach can sufficiently model any nonperiodic light curve. However, in the presence of periodicity, it is preferable to use a Seasonal (S) Autoregressive (AR) Integrated (I) Moving Average (MA) or SARIMA(p,d,q) × (P,D,Q)s model (Adhikari & Agrawal 2013), which can be expressed as

( 1 j = 1 p θ j L j ) ( 1 j = 1 P Θ j L sj ) Δ d Δ s D F ( t i ) = ( 1 + j = 1 p ϕ j L j ) ( 1 + j = 1 Q Φ j L sj ) ϵ ( t i ) , $$ \begin{aligned} \begin{aligned}&\implies \left( 1 - \sum _{j=1}^p \theta _j \mathcal{L} ^j \right) \left( 1 - \sum _{j=1}^{P} \Theta _j \mathcal{L} ^{sj} \right) \Delta ^d\Delta _s^D\mathcal{F} (t_i) \\&\quad = \left( 1 + \sum _{j=1}^p \phi _j \mathcal{L} ^j \right) \left( 1 + \sum _{j=1}^{Q} \Phi _j \mathcal{L} ^{sj} \right) \epsilon (t_i), \end{aligned} \end{aligned} $$(5)

where the second factor in both LHS and RHS is the seasonal term responsible for periodicity with a period s, ΔDℱ(ti) = (1−ℒs)Dℱ(ti) is the order of seasonal differencing, and P and Q are the orders of seasonal AR and MA, respectively.

Then the search for periodicity involves fitting the light curve using both ARIMA and SARIMA models and comparing their goodness-of-fit. We used the Akaike information criterion (AIC, Akaike 1974) for model comparisons. The AIC≡−2ln L+2k, where L is the likelihood of obtaining the data given the model and k is the number of free parameters in the model. AIC rewards a model for fitting the data better while penalizing it for using a larger number of parameters. We favor the model with the lowest AIC during model comparison. In the present case, if there is no periodicity, models with (p,d,q) × (0,0,0)s will have a lower AIC. We performed a grid search of AIC values in the parameter space

ϕ = { p , q [ 0 , 10 ] P , Q [ 0 , 6 ] d , D { 0 , 1 } s [ 0 , 10 ] days $$ \begin{aligned} \phi = {\left\{ \begin{array}{ll} p,q&\in [0,10]\\ P,Q&\in [0,6]\\ d,D&\in \{0,1\}\\ s&\in [0,10]\ \text{ days} \end{array}\right.} \end{aligned} $$(6)

and obtained the most likely model for the light curve. If a model with a periodic component better explains the light curve, the lowest AIC value will be in a model with nonzero P,  Q,  or D.

3.4. Significance estimation

Along with constructing the PSDs and identifying any dominant period, it is necessary to determine its significance. The significance of the PSD peak quantifies how likely it is to obtain a particular peak power due to random fluctuations, given the underlying model. This requires an assumption for the underlying model of the periodogram and, by extension, the underlying model of the light curve. The most basic model, ARIMA(0,0,0), assumes that each emission is independent, ℱ(ti) = ϵ(ti), which produces a PSD where the emitted power in different frequencies are independent (white-noise). We can model such a PSD using a constant in the frequency domain. With N different PSD frequencies, the probability (p) that the maximum power at any frequency crossing a threshold (z) is given by p(>z) ≈ N ez (Scargle 1982; Hong et al. 2018), quantifying the false alarm probability (FAP), or the significance of the peak. The lower the FAP of a period peak, the less likely it is to be caused by statistical fluctuation.

We empirically know that the underlying power spectrum of blazars (like most autoregressive process) demonstrates a red-noise periodogram model (e.g., Vaughan 2005), with more power being emitted at lower frequencies. So using an underlying ARIMA(0,0,0) or a white-noise model would give wrong estimates for the significance of the peaks. Ideally it makes sense to use the PSD model corresponding to the (S)ARIMA(p,d,q) model that best explains the given light curve. However, an analytical form for the periodogram is available only for ARIMA(1,0,0) model (or AR1 light curves, Robinson 1977), where the flux at a particular time depends only on the flux that preceded it. Percival & Walden (1993) models the ARIMA(1,0,0) PSD as

G rr ( f j ) = G 0 1 θ 2 1 2 θ cos ( π f i / f Nyq ) + θ 2 , $$ \begin{aligned} G_{rr}(f_j) = G_0\frac{1-\theta ^2}{1-2\theta \cos (\pi f_i/f_{\rm Nyq}) + \theta ^2}, \end{aligned} $$(7)

where fjs are the discrete frequencies up to the Nyquist frequency (fNyq), G0 is the average spectral amplitude, θ ≡ exp(Δt/τ) is the average autocorrelation coefficient and Δt is the average sampling interval. We obtain the characteristic autocorrelation timescale (τ) from Welch-overlapped-segment-averaging (WOSA, Welch 1967) of the LSP. We then estimated the significance of a peak by considering a χ2 distribution of periodogram values about the theoretical model. This procedure is performed by the computer code (REDFIT1, Schulz & Mudelsee 2002). The distribution is χ2 since the periodogram points are constructed from mod-squaring the real and imaginary parts which are assumed to be normally distributed.

Simpler light curves models, such as pure power-laws (P ∝ να1, Vaughan 2005) or smooth bending power-laws (Vaughan 2010) could reasonably approximate the underlying red-noise PSD. We observed (again using AIC) that the underlying model is closer to a smooth bending power-law than a simple power-law. So we used a smooth bending power-law to model the PSD for non-AR1 light curves. We used Monte Carlo (MC) techniques to simulate one thousand light curves with the same underlying PSD model and distribution of fluxes (PDF) as the original light curve (Emmanoulopoulos et al. 2013). From the LSP and WWZ analyses of the simulated light curves, we calculated the significance of the peak from the fraction of simulated light curves where the power crosses that of the original at the dominant frequency. We modeled the underlying PDF with a log-normal distribution.

4. Results

The three-hour binned γ-ray and R-band light curves from MJD 57710 to 57790 are presented in Fig. 1 along with the optimal block representation. We identify eight significant peaks during this period. From the DCF in Fig. 1, we also observe the γ and optical emissions to be significantly correlated at zero lag during this period with slightly lesser significant peaks at lags of ≈ + 7 days and ≈ − 14 days, giving us the first hint of periodicity. The WWZ map for this period shows strong periodicity in both the wavebands, spanning from MJD 57715 to 57770, with a dominant period of 7 . 6 0.17 + 0.31 $ 7.6^{+0.31}_{-0.17} $ in the γ-ray band and 7 . 6 0.25 + 0.36 $ 7.6^{+0.36}_{-0.25} $ in the optical R-band (Fig. 2). These values were obtained by fitting the most significant peak with a log-normal function. The same apparent quasi-periods were observed in the LSP of the two light curves in Fig. 2. To calculate the significance, the original periodogram was fit with a bending power-law model using a maximum likelihood method. It was confirmed (using AIC) that smooth bending power-law (AIC: 656) is a better fit for the periodograms than a regular power-law (AIC: 703). Using MC techniques from the thousands of simulated light curves, the significance of the dominant period is ∼4.0σ in both the wavebands. Considering the theoretical AR1 spectrum, dominant periods in both the wavebands show a significance of > 99%. The FAP for the dominant period in both the wavebands is < 10−5, although that considers an unrealistic underlying white-noise model (Fig. 2).

thumbnail Fig. 2.

WWZ maps of CTA 102 during the 2016–17 flare. (a) γ-ray WWZ map showing strong periodicity of ∼7.6 days. (b) Time averaged γ-ray WWZ map (black) and the LSP (red). The dominant period (∼7.6 days) shows > 5σ significance. (c) Optical R-band WWZ map showing similar periodicity. (d) Same as (b) for R-band emissions.

thumbnail Fig. 3.

LSP along with the theoretical model and line of 99% significance as generated by REDFIT. (a) For γ-ray emissions. (b) For R-band emissions.

While modeling the light curve using seasonal ARIMA, the model with the lowest AIC is SARIMA(4, 0, 0)×(0, 0, 1)7.625 d implying a periodic model better explains the light curve (Fig. 4). The SARIMA model was generated including an overall quadratic trend in the light curve during this period. The AIC values for the ARIMA models ignoring any periodic effect are given in Fig. 4a, while Fig. 4b gives the AIC values including seasonal effects for (p = 4, d = 0, q = 0) with a period of 7.625 days, where (P = 0, D = 0, Q = 1) gives the global minimum AIC. Figure 4c gives the AIC values considering different periods. We find the best AIC is for a period of 7.625 days which is consistent with our previous observations. The best model for the light curve during this period is then

F ( t i ) = j = 1 4 θ j F ( t i j ) + ϵ ( t i ) + ϵ ( t i 7.625 d ) . $$ \begin{aligned} \mathcal{F} (t_{i}) = \sum _{j=1}^4 \theta _j\ \mathcal{F} (t_{i-j}) + \epsilon (t_i) + \epsilon (t_{i-7.625\text{ d}}). \end{aligned} $$(8)

thumbnail Fig. 4.

AIC map for the SARIMA modeling of the light curve. (a) AIC map for ARIMA(p, d, q) models. Though the minimum occurs at ARIMA(3, 0, 2) (the box marked with magenta square), the global minimum, while considering SARIMA comes at SARIMA(4, 0, 0)×(0, 0, 1) (the box marked with red square). (b) AIC map for SARIMA(4, 0, 0)×(P, 0, Q). We observe the global minimum at SARIMA(4, 0, 0)×(0, 0, 1) considering a period (season of 7.625 days). The AIC values are provided in the boxes. In the uncolored region, either the fit did not converge or the AIC is greater than a set threshold. This is done to improve color resolution. (c) AIC values for different periods. We see a global minimum at 7.625 days.

It is best for ARIMA modeling if the data is evenly sampled; however, it can allow for a moderate amount of unsampled data (Feigelson et al. 2018). In the present case, we use ARIMA only on the γ-ray light curve which has nearly perfectly (98.4%) even sampling. We also tried linearly interpolating the for the small amount of unknown data and that did not change our conclusion.

We repeated the γ-ray analysis with bin sizes of 12 hours and 1 day to check any role of binning in the observed results. Both of these WWZ analyses demonstrated exactly the same periodicity with a significance of > 4σ.

5. Discussion

The results of the previous section (Sect. 4) demonstrate a periodic nature in the outburst of CTA 102 during 2016-end and 2017-beginning. The dominant period of ∼7.6 days is highly significant (> 4σ) as is evident from the numerous statistical tests performed. It is also to be noted that the γ-ray and R-band emissions are correlated and the dominant period is the same in both the wavebands. Also, the periodicity persists throughout the flaring period, with five QPO peaks that coincide with the significant oberved peaks (Fig. 1). Given that simultaneous activity was observed across this broad range of the accessible EM spectrum, it is extremely likely that the variability originated from the jet. Since the QPO was observed for almost the entire flaring period, understanding the origin of QPO could shed light on the origin of the flare itself.

Due to the transient nature of the quasi-periodic oscillation, we can disregard models that produce persistent QPOs such as a binary black hole system (Valtonen et al. 2008; Villforth et al. 2010) and persistent jet precession (Romero et al. 2000; Rieger 2004; Liska et al. 2018). Accretion-powered sources are known for exhibiting transient QPOs which are normally believed to be related to the accretion (Rodriguez et al. 2002). One possible model (for QPOs with tens of days period) involve hotspots rotating in the innermost stable circular orbit (ISCO) orbit of the SMBH (e.g., Zhang & Bao 1991; Gupta et al. 2009). Here the optical periodic flux modulation comes from the circular motion of the hotspot which modulates the seed photon field for the external Compton (EC) interaction with the particles in the jet, thereby producing modulation in the γ-ray flux (Gupta et al. 2009, 2019). This model is not favored, firstly since blazar emissions are mostly jet dominated in the optical band as well. Also, according to that model, the γ-ray emission is Doppler boosted and thus should not have the same dominant quasi-period as the optical.

The observed blazar emission, however, is dominated by the jet, which is ultimately powered by the accretion. It is thereby prudent to focus on the mechanisms that originate in the jet itself. QPOs related to jet are expected to be observed across the EM bands (e.g., Sarkar et al. 2020). On the contrary, QPOs due to binary interaction can be broadband like the jet or can be limited to only a few of the EM bands, depending on the interactions between the two and the mechanism responsible for the emission (e.g., Dey et al. 2018).

One possible model for this jet based periodicity is magnetic reconnection events in magnetic islands inside the jet (Huang et al. 2013). Multiple magnetic islands (X-point reconnection) situated roughly equidistantly will have their emissions delayed by a roughly constant time, giving rise to equidistantly spaced peaks in the observed flux, thereby mimicking rapid QPOs in standard BH systems. Shukla et al. (2018) observed extremely short timescale variability (∼5 min) in the source, which is smaller than the light travel time across the central SMBH (∼70 min); hence it was attributed to the magnetic reconnection events in the magnetic islands inside the jet.

Another, purely geometrical origin of the QPO, could be a region of enhanced emission (or blob) moving helically inside the jet (Mohan & Mangalam 2015; Sobacchi et al. 2017). A region of high particle density moving in the jet can enhance the jet luminosity. Due to such helical motion, the viewing angle (θobs) of the blob to line of sight varies as (Sobacchi et al. 2017; Zhou et al. 2018)

cos θ obs ( t ) = sin ϕ sin ψ cos ( 2 π t / P obs ) + cos ϕ cos ψ , $$ \begin{aligned} \cos \theta _{\rm obs}(t) = \sin \phi \ \sin \psi \cos (2\pi t/P_{\rm obs}) + \cos \phi \ \cos \psi , \end{aligned} $$(9)

where Pobs is the observed period, ψ is the viewing angle to the jet axis from our line of sight, and ϕ is the pitch angle of the helical motion. This variation in the viewing angle produces the periodic modulation in the Doppler factor (δ(t) = [Γ(1 − β cos θobs(t))]−1), which gives rise to the periodic modulation of the observed flux via F ν (t)= δ 3+α (t) F ν $ \mathcal{F}_{\nu}(t) = \delta^{3+\alpha}(t)\ \mathcal{F}^\prime_{\nu^\prime} $. Here F ν $ \mathcal{F}^\prime_{\nu^\prime} $ is the rest frame emission and α is the spectral index. This model can explain the transient nature of the periodicity. The QPO begins as the blob is injected near the base of the jet and it persists until the blob dissipates. Taking ϕ ≈ 2° (Zhou et al. 2018), ψ ≈ 3.7°, with the bulk Lorentz factor, Γ = 15.5 (Hovatta et al. 2009) and Pobs ≈ 7.6 days, the period in the rest frame of the blob would be

P = P obs ( 1 β cos ϕ cos ψ ) 4.3 years . $$ \begin{aligned} P=\frac{P_{\rm obs}}{(1-\beta \cos \phi \cos \psi )}\approx 4.3\ \text{ years}. \end{aligned} $$(10)

The total distance traveled by the blob during one period is D = cβP cos ϕ ≈ 1.31 pc. The projected distance traveled by the blob during the QPO observation is given by DP = 8D sin ψ ≈ 0.68 pc. Even though it did not occur in the same epoch, the short variability timescale ∼ 5 min (Shukla et al. 2018) can be explained in this scenario by considering a blob that is smaller than the central SMBH and has a very high γ. This would result in very high Doppler factor (δ) for very small θobs and the timescales of any intrinsic fluctuations (Δt) in the rest frame of the blob will then be compressed to (Δt/δ) in the observed frame, producing short variability timescales. The emission from a blob moving helically inside a straight jet will give rise to periodic modulations in the observed fluxes with nearly constant amplitude in each period. However, we do not observe this trend as the amplitudes are not constant throughout the periods. Much of this change in peak fluxes can possibly be explained by allowing for a curvature in the jet (Raiteri et al. 2017; Sarkar et al. 2020).

Of all the models discussed above, the most likely origin of the QPO, and by extension, the 2016–2017 flare, in CTA 102 is a sudden injection of a blob into the jet which then traverses in a helical motion. This conclusion is supported by an earlier claim of observations of a helical jet structure in the source by Li et al. (2018) in a nearby epoch where they use a helical jet to model radio variability and VLBA component trajectories. Also, Raiteri et al. (2017) explained the variability of the source using a curved jet which is essential to explain the modulation in the peak flux in different periods. We also observe a change in the optical polarization angle from 0° to 180° during the QPO period (Fig. 1c), which is an indication of helical motions in the jet. However, one difference is that Li et al. (2018) determined a theoretical period of ∼2 years which is much longer than the one observed, but can be explained by considering a higher bulk Lorentz factor. It must be kept in mind that blazar emission mechanisms are invariably complex, and processes like pulsational accretion flow instabilities often approximate periodic behaviors in the light curve (McKinney et al. 2012).

6. Conclusions

We report the detection of a significant multiwaveband (γ-ray and optical) QPO with week-like period in the flux of the blazar CTA 102 during its 2016–2017 outbursts. Claims for the detections of multiwaveband QPOs are rare in blazar light curves (to name a few Ackermann et al. 2015; Sarkar et al. 2020) as are putative QPOs with week- to month- like periods (Zhou et al. 2018; Sarkar et al. 2020). We stress that several independent techniques were employed to quantify the periodicity in the light curve and all of them produced consistent results of around ∼7.6d. Nonetheless, it is necessary to keep in mind that the significance of a QPO claim can depend on the parameters (bin sizes) and thresholds (test statistics or signal to noise ratio in each bin) adopted during the analysis. However, we found that the present γ-ray QPO detection is quite resilient to the parameters considered (as revealed from reanalyses with larger bin sizes). From our comparison of the data with models we conclude that the most likely origin of the flare (and the associated QPO) is a blob moving helically inside the relativistic jet.


Acknowledgments

We thank the anonymous referee for suggestions that improved the presentation of our results. This research has used data, software, and web tools of High Energy Astrophysics Science Archive Research Center (HEASARC), maintained by NASA’s Gooddard Space Flight Center. The work is also partly based on data collected by the WEBT collaboration and stored in the WEBT archive at the Osservatorio Astrofisico di Torino – INAF (http://www.oato.inaf.it/blazars/webt/); for questions regarding their availability, please contact the WEBT President Massimo Villata (massimo.villata@inaf.it). AS and VRC acknowledge support of the Department of Atomic Energy, Government of India, under project no. 12-R&D-TFR-5.02-0200. PK acknowledges an ARIES Aryabhatta Postdoctoral Fellowship (AO/A-PDF/770).

References

  1. Abdo, A. A., Ackermann, M., Agudo, I., et al. 2010, ApJ, 716, 30 [NASA ADS] [CrossRef] [Google Scholar]
  2. Acero, F., Ackermann, M., Ajello, M., et al. 2015, ApJS, 218, 23 [Google Scholar]
  3. Ackermann, M., Ajello, M., Albert, A., et al. 2015, ApJ, 813, L41 [NASA ADS] [CrossRef] [Google Scholar]
  4. Adhikari, R., & Agrawal, R. K. 2013, ArXiv e-prints [arXiv:1302.6613] [Google Scholar]
  5. Akaike, H. 1974, IEEE Trans. Autom. Control, 19, 716 [NASA ADS] [CrossRef] [MathSciNet] [Google Scholar]
  6. Atwood, W. B., Abdo, A. A., Ackermann, M., et al. 2009, ApJ, 697, 1071 [NASA ADS] [CrossRef] [Google Scholar]
  7. Bhatta, G. 2017, ApJ, 847, 7 [NASA ADS] [CrossRef] [Google Scholar]
  8. Bhatta, G. 2019, MNRAS, 487, 3990 [NASA ADS] [CrossRef] [Google Scholar]
  9. Bhatta, G., & Dhital, N. 2020, ApJ, 891, 120 [CrossRef] [Google Scholar]
  10. Blom, J. J., Bloemen, H., Bennett, K., et al. 1995, A&A, 295, 330 [NASA ADS] [Google Scholar]
  11. Brinkmann, W., Maraschi, L., Treves, A., et al. 1994, A&A, 288, 433 [Google Scholar]
  12. Carrasco, L., Dultzin-Hacyan, D., & Cruz-Gonzalez, I. 1985, Nature, 314, 146 [CrossRef] [Google Scholar]
  13. Casadio, C., Gómez, J. L., Jorstad, S. G., et al. 2015, ApJ, 813, 51 [NASA ADS] [CrossRef] [Google Scholar]
  14. Casadio, C., Marscher, A. P., Jorstad, S. G., et al. 2019, A&A, 622, A158 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  15. Cutri, R. M., Skrutskie, M. F., van Dyk, S., et al. 2003, VizieR Online Data Catalog: II/246 [Google Scholar]
  16. D’Ammando, F., Raiteri, C. M., Villata, M., et al. 2019, MNRAS, 490, 5300 [CrossRef] [Google Scholar]
  17. Dey, L., Valtonen, M. J., Gopakumar, A., et al. 2018, ApJ, 866, 11 [CrossRef] [Google Scholar]
  18. Edelson, R. A., & Krolik, J. H. 1988, ApJ, 333, 646 [NASA ADS] [CrossRef] [Google Scholar]
  19. Edelson, R. A., Saken, J., Pike, G., et al. 1991, ApJ, 372, L9 [CrossRef] [Google Scholar]
  20. Edelson, R., Krolik, J., Madejski, G., et al. 1995, ApJ, 438, 120 [NASA ADS] [CrossRef] [Google Scholar]
  21. Emmanoulopoulos, D., McHardy, I. M., & Papadakis, I. E. 2013, MNRAS, 433, 907 [NASA ADS] [CrossRef] [Google Scholar]
  22. Espaillat, C., Bregman, J., Hughes, P., & Lloyd-Davies, E. 2008, ApJ, 679, 182 [NASA ADS] [CrossRef] [Google Scholar]
  23. Feigelson, E. D., Babu, G. J., & Caceres, G. A. 2018, Front. Phys., 6, 80 [CrossRef] [Google Scholar]
  24. Fermi-LAT Collaboration 2020, ApJS, 247, 33 [Google Scholar]
  25. Foster, G. 1996, AJ, 112, 1709 [NASA ADS] [CrossRef] [Google Scholar]
  26. Fromm, C. M., Ros, E., Perucho, M., et al. 2013, A&A, 557, A105 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  27. Gasparyan, S., Sahakyan, N., Baghmanyan, V., & Zargaryan, D. 2018, ApJ, 863, 114 [NASA ADS] [CrossRef] [Google Scholar]
  28. Gaur, H., Gupta, A. C., Strigachev, A., et al. 2012, MNRAS, 425, 3002 [NASA ADS] [CrossRef] [Google Scholar]
  29. Gaur, H., Gupta, A. C., Bachev, R., et al. 2015, A&A, 582, A103 [CrossRef] [EDP Sciences] [Google Scholar]
  30. Grossmann, A., & Morlet, J. 1984, SIAM J. Math. Anal., 15, 723 [CrossRef] [MathSciNet] [Google Scholar]
  31. Gupta, A. C. 2014, JApA, 35, 307 [Google Scholar]
  32. Gupta, A. 2018, Galaxies, 6, 1 [CrossRef] [Google Scholar]
  33. Gupta, A. C., Fan, J. H., Bai, J. M., & Wagner, S. J. 2008, AJ, 135, 1384 [NASA ADS] [CrossRef] [Google Scholar]
  34. Gupta, A. C., Srivastava, A. K., & Wiita, P. J. 2009, ApJ, 690, 216 [NASA ADS] [CrossRef] [Google Scholar]
  35. Gupta, A. C., Tripathi, A., Wiita, P. J., et al. 2019, MNRAS, 484, 5785 [NASA ADS] [Google Scholar]
  36. Harris, D. E., & Roberts, J. A. 1960, PASP, 72, 237 [NASA ADS] [CrossRef] [Google Scholar]
  37. Hartman, R. C., Bertsch, D. L., Bloom, S. D., et al. 1999, ApJS, 123, 79 [NASA ADS] [CrossRef] [Google Scholar]
  38. Heidt, J., & Wagner, S. J. 1996, A&A, 305, 42 [NASA ADS] [Google Scholar]
  39. Hong, S., Xiong, D., & Bai, J. 2018, AJ, 155, 31 [NASA ADS] [CrossRef] [Google Scholar]
  40. Hovatta, T., Valtaoja, E., Tornikoski, M., & Lähteenmäki, A. 2009, A&A, 494, 527 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  41. Huang, C.-Y., Wang, D.-X., Wang, J.-Z., & Wang, Z.-Y. 2013, Res. Astron. Astrophys., 13, 705 [CrossRef] [Google Scholar]
  42. Kellermann, K. I., & Pauliny-Toth, I. I. K. 1968, ARA&A, 6, 417 [NASA ADS] [CrossRef] [Google Scholar]
  43. King, O. G., Hovatta, T., Max-Moerbeck, W., et al. 2013, MNRAS, 436, L114 [NASA ADS] [CrossRef] [Google Scholar]
  44. Kushwaha, P. 2020, Galaxies, 8, 15 [CrossRef] [Google Scholar]
  45. Kushwaha, P., Singh, K. P., & Sahayanathan, S. 2014, ApJ, 796, 61 [CrossRef] [Google Scholar]
  46. Kushwaha, P., Chandra, S., Misra, R., et al. 2016, ApJ, 822, L13 [NASA ADS] [CrossRef] [Google Scholar]
  47. Kushwaha, P., Sinha, A., Misra, R., Singh, K. P., & de Gouveia Dal Pino, E. M. 2017, ApJ, 849, 138 [CrossRef] [Google Scholar]
  48. Lachowicz, P., Gupta, A. C., Gaur, H., & Wiita, P. J. 2009, A&A, 506, L17 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  49. Larionov, V. M., Villata, M., Raiteri, C. M., et al. 2016, MNRAS, 461, 3047 [NASA ADS] [CrossRef] [Google Scholar]
  50. Li, X., Mohan, P., An, T., et al. 2018, ApJ, 854, 17 [NASA ADS] [CrossRef] [Google Scholar]
  51. Liska, M., Hesp, C., Tchekhovskoy, A., et al. 2018, MNRAS, 474, L81 [NASA ADS] [CrossRef] [Google Scholar]
  52. Lomb, N. R. 1976, Ap&SS, 39, 447 [Google Scholar]
  53. McKinney, J. C., Tchekhovskoy, A., & Bland Ford, R. D. 2012, MNRAS, 423, 3083 [NASA ADS] [CrossRef] [Google Scholar]
  54. Meyer, M., Scargle, J. D., & Blandford, R. D. 2019, ApJ, 877, 39 [NASA ADS] [CrossRef] [Google Scholar]
  55. Miller, H. R., Carini, M. T., & Goodrich, B. D. 1989, Nature, 337, 627 [NASA ADS] [CrossRef] [Google Scholar]
  56. Mohan, P., & Mangalam, A. 2015, ApJ, 805, 91 [NASA ADS] [CrossRef] [Google Scholar]
  57. Moore, R. L., & Stockman, H. S. 1981, ApJ, 243, 60 [NASA ADS] [CrossRef] [Google Scholar]
  58. Osterman Meyer, A., Miller, H. R., Marshall, K., et al. 2009, AJ, 138, 1902 [NASA ADS] [CrossRef] [Google Scholar]
  59. Percival, D. B., & Walden, A. T. 1993, Spectral Analysis for Physical Applications (Cambridge University Press) [CrossRef] [Google Scholar]
  60. Pica, A. J., Smith, A. G., Webb, J. R., et al. 1988, AJ, 96, 1215 [NASA ADS] [CrossRef] [Google Scholar]
  61. Quirrenbach, A., Witzel, A., Wagner, S., et al. 1991, ApJ, 372, L71 [NASA ADS] [CrossRef] [Google Scholar]
  62. Raiteri, C. M., Villata, M., Tosti, G., et al. 2003, A&A, 402, 151 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  63. Raiteri, C. M., Villata, M., Acosta-Pulido, J. A., et al. 2017, Nature, 552, 374 [NASA ADS] [CrossRef] [Google Scholar]
  64. Rieger, F. M. 2004, ApJ, 615, L5 [NASA ADS] [CrossRef] [Google Scholar]
  65. Robinson, P. 1977, Stochastic Processes Appl., 6, 9 [CrossRef] [Google Scholar]
  66. Rodriguez, J., Varnière, P., Tagger, M., & Durouchoux, P. 2002, A&A, 387, 487 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  67. Romero, G. E., Chajet, L., Abraham, Z., & Fan, J. H. 2000, A&A, 360, 57 [NASA ADS] [Google Scholar]
  68. Sandage, A., & Wyndham, J. D. 1965, ApJ, 141, 328 [NASA ADS] [CrossRef] [Google Scholar]
  69. Sandrinelli, A., Covino, S., & Treves, A. 2014, ApJ, 793, L1 [NASA ADS] [CrossRef] [Google Scholar]
  70. Sandrinelli, A., Covino, S., & Treves, A. 2016, ApJ, 820, 20 [NASA ADS] [CrossRef] [Google Scholar]
  71. Sandrinelli, A., Covino, S., Treves, A., et al. 2017, A&A, 600, A132 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  72. Sarkar, A., Gupta, A. C., Chitnis, V. R., & Wiita, P. J. 2020, ArXiv e-prints [arXiv:2003.01911] [Google Scholar]
  73. Scargle, J. D. 1981, ApJS, 45, 1 [NASA ADS] [CrossRef] [Google Scholar]
  74. Scargle, J. D. 1982, ApJ, 263, 835 [Google Scholar]
  75. Scargle, J. D., Norris, J. P., Jackson, B., & Chiang, J. 2013, ApJ, 764, 167 [NASA ADS] [CrossRef] [Google Scholar]
  76. Scaringi, S., Maccarone, T. J., Kording, E., et al. 2015, Sci. Adv., 1, e1500686 [NASA ADS] [CrossRef] [Google Scholar]
  77. Schmidt, M. 1965, ApJ, 141, 1295 [NASA ADS] [CrossRef] [Google Scholar]
  78. Schulz, M., & Mudelsee, M. 2002, Comput. Geosci., 28, 421 [NASA ADS] [CrossRef] [Google Scholar]
  79. Sembay, S., Warwick, R. S., Urry, C. M., et al. 1993, ApJ, 404, 112 [NASA ADS] [CrossRef] [Google Scholar]
  80. Shukla, A., Mannheim, K., Patel, S. R., et al. 2018, ApJ, 854, L26 [NASA ADS] [CrossRef] [Google Scholar]
  81. Sillanpaa, A., Takalo, L. O., Pursimo, T., et al. 1996, A&A, 305, L17 [Google Scholar]
  82. Sobacchi, E., Sormani, M. C., & Stamerra, A. 2017, MNRAS, 465, 161 [NASA ADS] [CrossRef] [Google Scholar]
  83. Stein, W. A., Odell, S. L., & Strittmatter, P. A. 1976, ARA&A, 14, 173 [NASA ADS] [CrossRef] [Google Scholar]
  84. Urry, C. M., & Padovani, P. 1995, PASP, 107, 803 [NASA ADS] [CrossRef] [Google Scholar]
  85. Urry, C. M., Maraschi, L., Edelson, R., et al. 1993, ApJ, 411, 614 [NASA ADS] [CrossRef] [Google Scholar]
  86. Valtaoja, E., Lehto, H., Teerikorpi, P., et al. 1985, Nature, 314, 148 [CrossRef] [Google Scholar]
  87. Valtonen, M. J., Lehto, H. J., Nilsson, K., et al. 2008, Nature, 452, 851 [NASA ADS] [CrossRef] [PubMed] [Google Scholar]
  88. Vaughan, S. 2005, A&A, 431, 391 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  89. Vaughan, S. 2010, MNRAS, 402, 307 [NASA ADS] [CrossRef] [Google Scholar]
  90. Véron-Cetty, M. P., & Véron, P. 2010, A&A, 518, A10 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  91. Villforth, C., Nilsson, K., Heidt, J., et al. 2010, MNRAS, 402, 2087 [Google Scholar]
  92. Wagner, S. J., & Witzel, A. 1995, ARA&A, 33, 163 [NASA ADS] [CrossRef] [Google Scholar]
  93. Welch, P. 1967, IEEE Trans. Audio Electroacoust., 15, 70 [Google Scholar]
  94. Welsh, W. F. 1999, PASP, 111, 1347 [Google Scholar]
  95. White, R. J., & Peterson, B. M. 1994, PASP, 106, 879 [NASA ADS] [CrossRef] [Google Scholar]
  96. Zacharias, M., Böttcher, M., Jankowsky, F., et al. 2019, ApJ, 871, 19 [NASA ADS] [CrossRef] [Google Scholar]
  97. Zhang, X. H., & Bao, G. 1991, A&A, 246, 21 [Google Scholar]
  98. Zhang, P.-F., Yan, D.-H., Zhou, J.-N., et al. 2017, ApJ, 845, 82 [NASA ADS] [CrossRef] [Google Scholar]
  99. Zhou, J., Wang, Z., Chen, L., et al. 2018, Nat. Commun., 9, 4599 [CrossRef] [Google Scholar]

All Figures

thumbnail Fig. 1.

CTA 102 light curves during the 2016–17 flare. (a) Fermi-LAT (0.1−300 GeV) light curve. The black line gives the optimal block representation of the light curve. A peak is considered to be significant if the block height is greater than a set threshold. (b) Optical R band light curve and its optimal block representation. The QPO peaks are denoted by red arrows. (c) γ-ray and R-band emission DCF. (d) Optical polarization angle during the detected QPO.

In the text
thumbnail Fig. 2.

WWZ maps of CTA 102 during the 2016–17 flare. (a) γ-ray WWZ map showing strong periodicity of ∼7.6 days. (b) Time averaged γ-ray WWZ map (black) and the LSP (red). The dominant period (∼7.6 days) shows > 5σ significance. (c) Optical R-band WWZ map showing similar periodicity. (d) Same as (b) for R-band emissions.

In the text
thumbnail Fig. 3.

LSP along with the theoretical model and line of 99% significance as generated by REDFIT. (a) For γ-ray emissions. (b) For R-band emissions.

In the text
thumbnail Fig. 4.

AIC map for the SARIMA modeling of the light curve. (a) AIC map for ARIMA(p, d, q) models. Though the minimum occurs at ARIMA(3, 0, 2) (the box marked with magenta square), the global minimum, while considering SARIMA comes at SARIMA(4, 0, 0)×(0, 0, 1) (the box marked with red square). (b) AIC map for SARIMA(4, 0, 0)×(P, 0, Q). We observe the global minimum at SARIMA(4, 0, 0)×(0, 0, 1) considering a period (season of 7.625 days). The AIC values are provided in the boxes. In the uncolored region, either the fit did not converge or the AIC is greater than a set threshold. This is done to improve color resolution. (c) AIC values for different periods. We see a global minimum at 7.625 days.

In the text

Current usage metrics show cumulative count of Article Views (full-text article views including HTML views, PDF and ePub downloads, according to the available data) and Abstracts Views on Vision4Press platform.

Data correspond to usage on the plateform after 2015. The current usage metrics is available 48-96 hours after online publication and is updated daily on week days.

Initial download of the metrics may take a while.