Intra-Night Variability of OJ 287 with Long-Term Multiband Optical Monitoring

Wei Zeng 1,2 ID , Qing-Jiang Zhao 3, Ze-Jun Jiang 2,4, Zhi-Hui Kong 4, Zhen Liu 1,2, Dong-Dong Wang 1,2, Xiong-Fei Geng 1,2, Shen-Bang Yang 1,2 and Ben-Zhong Dai 1,2,* ID 1 Department of Physics, Yunnan University, Kunming 650091, China; zw.z@me.com (W.Z.); zheng_liu1@hotmail.com (Z.L.); wdd92@outlook.com (D.-D.W.); gxfqjsf@163.com (X.-F.G.); shb_yang@hotmail.com (S.-B.Y.) 2 Key Laboratory of Astroparticle Physics of Yunnan Province, Yunnan University, Kunming 650091, China 3 Department of Physics and Technology, Kunming College, Kunming 650214, China; zhaoqj213@163.com 4 Department of Astronomy, Yunnan University, Kunming 650091, China; zjjiang@ynu.edu.cn (Z.-J.J.); kongzh0103@163.com (Z.-H.K.) * Correspondence: bzhdai@ynu.edu.cn


Introduction
BL Lacertae (BL Lac) objects, which have either very weak or no emission lines [1], and flatspectrum radio quasars (FSRQs) with strong emission lines [2,3] form a subclass of radio-loud active galactic nuclei (AGNs) known as blazars. Blazars are characterised by non-thermal emission, and strong and rapid flux variability across the entire electromagnetic spectrum, from radio-to gamma-rays. The emission is normally attributed to the relativistic jet oriented at a small angle from the line of sight [4]. Blazar flux variability timescales extending from a few minutes to years and even decades can be broadly divided into three classes: a large variation over hours to days is often known as intra-day or -night variability (IDV) or micro-variability [5,6]; variation on timescales of days to weeks, or even a few months, is considered short-term variability (STV); meanwhile, long-term variability (LTV) can have timescales from several months to years [7][8][9][10].
OJ 287 (α = 08:54:48.9, δ = +20:06:30.6, J2000) is a blazar at redshift z = 0.306 [11]. Sillanpää et al. [12] pointed out for the first time that there is a double-peak structure in the cyclic optical outbursts of OJ 287 by using the optical V-band observations starting from 1890. The curve exhibited periodic outbursts at intervals of ∼12 years. A binary pair of super-massive black holes system was used to explain this quasi-periodic light curve. Sillanpää et al. [13] (1996a [13], 1996b [14]) reported the next double-peak outburst occurrence in 1994-1995, which occurred almost exactly at the predicted times. Then a double-peaked outburst was also seen during the next recurrence in 2005-2008 [15]. These predicted recurrences are usually interpreted as OJ 287 housing a binary pair of black holes system with a period of ∼12 years.
In the case of OJ 287, outbursts occurring roughly every 12 years are almost certainly produced, and it was predicted that OJ 287 should show a major outburst in 2015-2016 [16]. In our work, along with the confirmation of this predicted outburst during 2015-2016 reported by [10,17], we study the properties of variability and spectral variation in a short timescales.

Observations and Data Reduction
Our photometric observations of blazar OJ 287 were performed using the 1.02 m (YAO1.02) optical telescopes at the Yunnan Astronomical Observatory. Since 2009, this telescope has been equipped with an Andor DW436 CCD camera (2048 × 2048 pixels) at the Cassegrain focus ( f = 13.3 m). The readout noise and gain are 6.33 electrons and 2.0 electrons/ADU in a 2 µs readout time model, or 2.29 electrons and 1.4 electrons/ADU in a 16 µs readout time model. The FOV of the CCD frame is 7.3 × 7.3 arcmin 2 with a pixel scale of 0.21 arcsec/pixel. We used standard Johnson broad-band filters with effective wavelength midpoints of V = 551 nm, R = 658 nm, and I = 806 nm.
We performed optical observations in the V-, R-, and Ibands in the cyclic mode. The exposure times from 50 to 300 s were chosen according to seeing, the filter and the telescope. The time resolutions (the time interval between two adjacent data points in the same band) were less the 15 min and were approximately 10 min in most cases. Therefore, these data in the cyclic mode were considered as quasi-simultaneous measurements, which were explored for analysing inter-band correlation and colour index. During our observation campaign (6 March 2010 to 3 April 2016), we observed for a total of 34 nights, to obtain 2255 CCD frames (shown in Table 1). Table 1. A short summary of the observations in each band performed from 2010 and 2016. In total we observed for 34 nights and obtained 2255 CCD frames. 2010  10  183  196  197  2012  4  133  136  136  2013  7  142  153  156  2014  5  42  46  55  2015  2  19  25  39  2016  6  199  201  197 All images were reduced with standard image reduction and analysis facility (IRAF) 1 procedures after bias and flat-field corrections. Aperture photometry was performed on the source and comparison stars with APPHOT. Photometry of the source and comparison stars was performed with the same aperture, which was determined by the full width at half-maximum (FWHM) of the comparison stars and was the same for each observation. We compared the photometric results and found that an aperture radius of about 1.5 FWHM almost always provided the best S/N ratio, and thus we concentrated on this aperture for our analysis. Each night, we performed aperture photometry with different aperture radii.

Year Nights V(N) R(N) I(N)
For each CCD frame, the instrumental magnitudes of OJ 287 and the two comparison stars (listed in Table 2) were extracted. The brightness (magnitude) of OJ 287 was calculated as the average of that derived with respect to comparison stars 4 and 10, and the corresponding standard deviation was treated as the error of OJ 287 within each CCD frame. The deviation of the average of the differential instrumental magnitude of comparison stars 4 and 10, delta (star 4-star 10), was used to verify the 1 IRAF is distributed by the National Optical Astronomy Observatories, which are operated by the Association of Universities for Research in Astronomy, Inc., under cooperative agreement with the National Science Foundation. stability of the comparison stars, which is also taken as the accuracy of the observations. Figure 1 shows the deviation of comparison stars 4 and 10 for the entire observational campaign. The deviations in bands V-, Rand I-bands are mainly distributed in ±0.025 mag.

Variability
The V, R and I light curves from 6 March 2010 to 3 April 2016 are shown in Figure 2. The average magnitudes in each band are V = 14.728 ± 0.352, R = 14.351 ± 0.328 and I = 13.763 ± 0.309. The amplitude variability (maximum minus minimum) in each filter are respectively ∆V = 1. m 335, ∆R = 1. m 237 and ∆I = 1. m 594.
Presently, a number of statistical tests, such as the C-test, the F-test, the χ 2 -test and one-way analysis of variance (ANOVA), have been proposed to assess quantitatively whether there are IDVs [6,9,[19][20][21][22][23]. Here, we apply the F-test and χ 2 -test to cross-check the intra-day light curves. We consider variability when the light curve satisfied both criteria described below.
The F-test is regarded as a powerful distribution statistic to check for the presence of variability, as introduced by [21]. When comparing two sample variances, the F-statistic values are calculated as: where S 2 BL−StarA , S 2 BL−StarB and S 2 StarA−StarB are the variances of the differential instrumental magnitudes of the blazar and comparison star A, the blazar and comparison star B, and comparison stars A and B, respectively. The F-statistic is compared with the F (α) ν BL ,ν c critical value. The number of degrees of freedom for each sample, ν BL and ν c , are the number of measurements minus 1 (N − 1), and α is the significance level set for the test. In this paper, the F-test was performed at two significance levels: 0.99 and 0.999. If either F 1 or F 2 was larger than the critical value, the null hypothesis (no variability) was discarded. The χ 2 -test was also used: where m i are the individual magnitudes, σ m i are their errors, and m denotes the mean value of magnitude for an observation night. This statistic is also compared against the critical value χ 2 α,ν obtained from the χ 2 probability function, where ν is the degrees of freedom and α is the significance level. χ 2 > χ 2 α,ν provides evidence of variability. OJ 287 was said to be variable (V) if the statistic from both tests was satisfied at the 0.999 level, while it was marked as "non-variable" (NV) if none of these criteria were met at the 0.99 level. For neither, it was marked as "probably variable" (PV). The results of the Fand χ 2 -test and the IDV observations in the V-, Rand I-bands are presented in Table 3. As shown in Table 3, 17 light curves were marked as variable among 100 light curves (18 light curves met the criterion with the F-test and 57 light curves with the χ 2 -test). As an illustration, the significant variability on two nights are presented, 6-8 January in 2016, which were marked as variable in the V-, Rand I-bands (shown in Figure 3). On January 6, 2016, the magnitudes showed a monotonically increasing trend. The magnitude changes of the V-, Rand I-bands were ∆V = 0.    For each IDV, the variability amplitude (Amp) could be calculated by [24]: where A max and A min are the maximum and minimum magnitude, respectively, of the light curve for the night being considered, and σ is the corresponding standard deviation. Our observational data were obtained with different exposures and time intervals across 6 years. To avoid uncertainties caused by different observational set-ups, the fractional variability amplitude (F var ) was calculated to show the intrinsic variability amplitude of the source by removing the effects of measurement noise; it is defined as follows [25]: Here S 2 denotes the total variance of the light curve: where σ 2 err is the mean error squared, and m is the mean value of the magnitude. The error on F var is Our results of the IDV analysis are shown in Figure 4 and Table 3. Both the statistical distributions of two variability amplitudes show that the IDV amplitudes of OJ 287 during the whole observational campaign were small.   Table 3.
We extracted the variability for candidate light curves (marked as PV or V states) by an asymmetric flare template [26], which describes flares as exponentially rising and decaying by the function: where F 0 is the constant background flux level (could be a constant or linear model), F i is the flare normalisation, t i is the time of maximum flare, T r,i is the flux rising timescale, and T d,i is the flux decaying timescale. Most of candidate light curves showed a monotonically increasing or decreasing trend, which meant the timescale of the variability was larger than the length of the data. Only two V light curves and four PV light curves showed complete flares. The results of the light curves are shown in Figure 5 and Table 4.

Colour Index
It is common that brightness variations are often associated with changes in spectral shapes. Thus we also investigated the correlations between the colour indices and magnitude variations.
Here we consider whether variations in the V − R colour indices of OJ 287 changed with respect to variations in its brightness in the R-band on short-and long-term bases. These colour and magnitude plots of OJ 287 are displayed in Figures 6 and 7. We calculated the best linear fit as shown by the straight lines in Figures 6 and 7 for the colour index against magnitude, and the slope, intercept, linear Pearson correlation coefficients and corresponding null hypothesis probability values are also listed in Table 5. A positive Pearson's r coefficient between the colour index and apparent magnitude of the blazar implies a positive correlation, which means the source tends to be bluer when it is brighter (BWB) or redder when fainter (RWF), while a negative Pearson's r coefficient suggests the opposite correlation: redder when brighter (RWB) or bluer when fainter (BWF) behaviour. We found a weak positive correlation between the colour indices and the R-band magnitude on long timescales during these observations of OJ 287 (shown as Figure 6). Then we measured Pearson's r coefficient between the colour indices and the R-band magnitude on intra-night timescales, as plotted in Figure 7. Ten nights showed a significant negative correlation between the colour indices and R-band magnitude on intra-night timescales, with probability values below 0.001; other nights showed nearly achromatic (ACH) behaviour.    Table 5.

Discussion and Conclusions
Blazars, being dominated by a relativistic jet oriented at a small angle from the line of sight, display variability ranging from a few hundredths of a magnitude to more than several magnitudes. Their variability timescales range from minutes to years. Either factors intrinsic to the jet or due to changes in overall jet power could explain such variability. The shock-in-jet model can be used to explain most of variabilities in blazars [10,[27][28][29][30]. Any changes in the physical mechanisms of blazars, such as to the magnetic field, electron density or velocity, can trigger a shock that leads to flares when propagating along the relativistic jet. Turbulence behind the shock-in-jet can be a good way to explain the fast variability in the IDV of blazars [10,[31][32][33]. However, when blazars are in a low state, the emission from the accretion disk dominates the jet emission. Instability in the accretion disk, such as at "hot spots", can also produce IDVs [9,10,[34][35][36]. Perturbations on the accretion disks transfer into the relativistic jet and, being Doppler amplified, might also produce fluctuations [9,10,37,38].
Previous colour index studies suggest apparent dichotomy responses: BL Lacs become BWB and FSRQs become RWB [9,10,[39][40][41][42][43][44]. The BWB behavior in BL Lacs can be interpreted as electrons being accelerated to preferentially higher energies before radiatively cooling, and the RWB behavior in FSRQs is explained by the addition of redder, non-thermal jet emission to an already bluer, thermal disk component. These results are often preferentially observed during flaring, which causes the selection effect that the two extreme ends of the colour-magnitude correlation are most regularly reported. Ruan et al. [45] point out that spectral variability could be the result of hot spots in the accretion disk emission by using a sample of 604 variable quasars. Isler et al. [46] suggest colour-magnitude variability should be a continuum rather than a dichotomy. They use empirical data from 3C 279 to explain how blazars can smoothly evolve from a jet-quiescent, disk-dominated colour profile to an actively jet-dominated state and associated colour profile and back to a jet-quiescent state. Our results of colour-magnitude correlations on long and short timescales tend to be interpreted in the framework of a relatively bluer accretion disk and redder jet change in intensity. In addition, the jet intensity could change in a short timescales, while the disk intensity changed more slowly. Thus the intra-night colour behavior was alternatively BWF, RWB or ACH, which can be described by the variation in contributions of the non-thermal relativistic jet in OJ 287.
We carried out multi-band optical photometric monitoring of the blazar OJ 287 over 34 observation nights between 6 March 2010 and 3 April 2016. We searched for flux and colour-magnitude variation on IDV timescales. The IDV amplitudes of OJ 287 were small for the observed flux variabilities. The blazar was in a highly active and variable state during January 2016, showing relatively large amplitude flux variations: F var of 1.3−2.1%, and Amp of 5.8−9.0%. OJ 287 showed a weak positive correlation on the long timescale, while negative correlations were found on IDV timescales. Our present analysis cannot distinguish between the possible physical mechanisms in detail. Very dense and highly precise simultaneous multi-band observations are necessary. Thus OJ 287 should continue to be monitored whenever possible.