arXiv:1505.07467v2 [astro-ph.HE] 22 Mar 2016
MNRAS
000
,
1
–
16
(2015)
Preprint 10 August 2016
Compiled using MNRAS L
A
T
E
X style file v3.0
RoboPol: first season rotations of optical polarization pla
ne
in blazars
D. Blinov
1
,
7
⋆
, V. Pavlidou
1
,
2
, I. Papadakis
1
,
2
, S. Kiehlmann
3
, G. Panopoulou
1
,
I. Liodakis
1
, O. G. King
4
, E. Angelakis
3
, M. Balokovi ́c
4
, H. Das
6
, R. Feiler
5
,
L. Fuhrmann
3
, T. Hovatta
8
, P. Khodade
6
, A. Kus
5
, N. Kylafis
2
,
1
, A. Mahabal
4
,
I. Myserlis
3
, D. Modi
6
, B. Pazderska
5
, E. Pazderski
5
, I. Papamastorakis
1
,
2
,
T. J. Pearson
4
, C. Rajarshi
6
, A. Ramaprakash
6
, P. Reig
2
,
1
, A. C. S. Readhead
4
,
K. Tassis
1
,
2
, J. A. Zensus
3
1
Department of Physics and Institute for Plasma Physics, Uni
versity of Crete, GR-71003, Heraklion, Greece
2
Foundation for Research and Technology - Hellas, IESL, Vout
es, 7110 Heraklion, Greece
3
Max-Planck-Institut f
̈
ur Radioastronomie, Auf dem H
̈
ugel 69, 53121 Bonn, Germany
4
Cahill Center for Astronomy and Astrophysics, California I
nstitute of Technology, 1200 E California Blvd, MC 249-17,
Pasadena CA, 91125, USA
5
Toru ́n Centre for Astronomy, Nicolaus Copernicus Universi
ty, Faculty of Physics, Astronomy and Informatics,
Grudziadzka 5, 87-100 Toru ́n, Poland
6
Inter-University Centre for Astronomy and Astrophysics, P
ost Bag 4, Ganeshkhind, Pune - 411 007, India
7
Astronomical Institute, St. Petersburg State University,
Universitetsky pr. 28, Petrodvoretz, 198504 St. Petersbur
g, Russia
8
Aalto University Mets
̈
ahovi Radio Observatory, Mets
̈
ahovintie 114, FL-02540 Kylm
̈
al
̈
a, Finland
Accepted XXX. Received YYY; in original form ZZZ
ABSTRACT
We present first results on polarization swings in optical emission of b
lazars obtained
by RoboPol, a monitoring programme of an unbiased sample of gamma-
ray bright
blazars specially designed for effective detection of such events. A
possible connection
of polarization swing events with periods of high activity in gamma rays
is investigated
using the data set obtained during the first season of operation. I
t was found that
the brightest gamma-ray flares tend to be located closer in time to r
otation events,
which may be an indication of two separate mechanisms responsible fo
r the rotations.
Blazars with detected rotations during non-rotating periods have
significantly larger
amplitude and faster variations of polarization angle than blazars wit
hout rotations.
Our simulations show that the full set of observed rotations is not a
likely outcome
(probability
≤
1
.
5
×
10
−
2
) of a random walk of the polarization vector simulated by
a multicell model. Furthermore, it is highly unlikely (
∼
5
×
10
−
5
) that none of our
rotations is physically connected with an increase in gamma-ray activ
ity.
Key words:
polarization – galaxies: active – galaxies: jets – galaxies: nuclei
1 INTRODUCTION
Blazars are active galactic nuclei whose jets are oriented
close to our line of sight, so that we observe high relativis-
tic beaming of their non-thermal emission and large ampli-
tude variability at all wavelengths. The low-frequency emi
s-
sion is dominated by synchrotron radiation, and hence is
highly polarized. The exact polarization fraction and di-
rection depend on the structure of the magnetic field in
the emitting region, and on the number of emitting re-
⋆
E-mail: blinov@physics.uoc.gr
gions along the line of sight. The polarization direction
(in the simple case of a single dominant emission region)
traces (and is perpendicular to) the direction of the pro-
jected magnetic field on the plane of the sky. Already from
early optical observations, it has been known that polariza
-
tion parameters of blazars are variable on daily time-scale
s
(
Kinman et al. 1966
). In general, both flux density and po-
larization exhibit an erratic variability (
Angel & Stockman
1980
;
Uemura et al. 2010
;
Ikejiri et al. 2011
), which could
be interpreted as a random walk (
Moore et al. 1982
). How-
ever, in some cases the electric vector position angle (EVPA
)
of the polarized emission displays long, smooth and mono-
c
2015 The Authors
2
D. Blinov et al.
Table 1.
Selection criteria for the gamma-ray–loud and the control s
ample.
Property
Gamma-ray–loud sample
Control sample
2FGL
included
not included
2FGL F(E
>
100 MeV)
>
10
−
8
cm
−
2
s
−
1
−
2FGL source class
agu, bzb, or bzq
−
Galactic latitude
|
b
|
>
10
◦
−
Elevation (Elv) constraints
1
Elv
max
≥
40
◦
for at least 90 consecutive
nights in the window June – November
Elv
max
≥
40
◦
for at least 90 consecutive
nights in the window April – November
R
magnitude
2
≤
17
.
5
≤
17
.
5
CGRaBS/15 GHz OVRO monitoring
no constraints
included
OVRO 15 GHz mean flux density
no constraints
≥
0
.
060 Jy
OVRO 15 GHz intrinsic modulation in-
dex,
m
no constraints
≥
0
.
05
1
Refers to elevation during Skinakas dark hours
2
Average value between archival value and measured during pr
eliminary RoboPol observations (when applicable)
tonic rotations which have been observed in the optical sinc
e
the 1980s (
Kikuchi et al. 1988
). A number of mechanisms
have been proposed for the interpretation of such events, in
-
cluding: stochastic variations of turbulent magnetic field
s,
a shock travelling through a non-axisymmetric magnetic
field (
Konigl & Choudhuri 1985
), polarized flares in the ac-
cretion disc (
Sillanp
̈
a
̈
a et al. 1993
), two-component models
consisting of two independent sources of polarized emissio
n
(
Bjornsson 1982
), and jet bending (
Abdo et al. 2010a
).
Blazars represent the most common class of known
gamma-ray sources (
Nolan et al. 2012
;
Acero et al. 2015
).
Despite the recent progress in the field, many questions con-
cerning the high-energy emission produced by blazars are
still under debate. For instance, it is unclear where the
gamma-ray emitting site is located: within the broad-line
region (e.g.
Blandford & Levinson 1995
;
Poutanen & Stern
2010
) or well downstream in the jet (e.g.
Marscher et al.
2008
;
Agudo et al. 2011
).
Recent work showed that at least some large EVPA
swings can be associated with gamma-ray flares (e.g.
Abdo et al. 2010a
;
Larionov et al. 2013b
) and therefore can
possibly provide some insight on the physics of high-energy
activity. Although such events have triggered an increasin
g
interest in polarimetric monitoring of gamma-ray blazars,
efforts in this direction have been based on selected cases
comprising statistically biased samples. As a result, a sig
-
nificant amount of invaluable polarimetric data sets for a
large number of sources has been gathered. However, this set
cannot be used for statistically rigorous population studi
es
and, in particular, the investigation of a possible correla
tion
between gamma-ray flares and optical EVPA rotations. The
RoboPol programme (
King et al. 2014
;
Pavlidou et al. 2014
)
has been designed to provide a data set of rotation events in
an unbiased sample of blazars, appropriate for such studies
.
In this paper, we analyse EVPA rotations detected by
RoboPol during the first observing season between 2013 July
and November. After a brief description of observing and
reduction techniques in Sec.
2
, we estimate the frequency of
EVPA rotations in blazars and list their properties in Sec.
3
.
A Monte Carlo simulation is performed in Sec.
4
in order
to determine whether the EVPA rotations can be produced
by random walk processes. In Sec.
5
we study the possible
connection between an increased activity in the gamma-ray
band and EVPA swings. Our findings are summarized in
Sec.
6
.
2 OBSERVATIONS AND DATA REDUCTION
2.1 Our sample
A unique feature of the RoboPol programme is that it is
monitoring a sample which has been selected on the basis of
strict, bias-free and objective criteria (for detailed dis
cussion
on the sample construction, see
Pavlidou et al. 2014
). The
sample consists of three distinct groups.
(i) The main (“gamma-ray–loud”) sample is an unbi-
ased subset of a statistically complete flux-limited sample
of blazars from the second
Fermi
-LAT source catalogue
(
Nolan et al. 2012
). Specifically, we selected all the sources
in the 2FGL catalogue classified as BL Lacertae objects
(bzb), Flat Spectrum Radio Quasars (bzq), or active galaxy
of uncertain type (agu). Applying the selection criteria li
sted
in Table
1
, we constructed a gamma-ray flux-limited“parent
sample”. Application of the visibility constraints and fiel
d-
quality cuts resulted in an unbiased subsample of 83 sources
,
among which we randomly selected 62 sources.
(ii) A “control” sample of 15 “gamma-ray–quiet” sources.
It constitutes an unbiased subset of a statistically comple
te
sample of blazars. It has been drawn from the CGRaBS
catalogue (
Healey et al. 2008
) applying the selection criteria
listed in Table
1
.
(iii) 24 additional sources chosen on the basis of their
variability characteristics or their presence either in th
e F-
GAMMA programme sample or in TeV catalogues.
MNRAS
000
,
1
–
16
(2015)
RoboPol: EVPA rotations in blazars
3
Figure 1.
Season length and median cadence for the first season
data. Broken lines limit areas where rotations at rates of 10
, 15
and 20 deg d
−
1
can be detected (see Sec.
3.3
for details). Only
objects with ∆
t
-median
≤
20 were left for the more detailed view.
Although here we present the polarization swings detected
in all monitored sources during the first RoboPol observing
season, the statistical analysis in this paper is based only
on
sources from the group (i).
2.2 Optical observations
All photometric and polarimetric measurements were done
at the 1.3-m telescope of Skinakas observatory
1
using
RoboPol, a polarimeter specifically built for the project
(
King et al. 2014
). The RoboPol instrument contains a fixed
set of two Wollaston prisms and half-wave plates, which
splits each incident ray into four rays with polarization
plane rotated 45
◦
with respect to each other. Measuring
relative intensities in pairs of the rays for each object in
the 13 arcmin
×
13 arcmin field, we obtain the fractional
Stokes parameters
q
= (
I
1
−
I
2
)
/
(
I
1
+
I
2
) =
Q/I
and
u
= (
I
3
−
I
4
)
/
(
I
3
+
I
4
) =
U/I
. Stokes parameter
I
is cal-
culated as a sum of intensities of all four spots. Since the
polarization parameters are measured simultaneously, we
avoid unmeasurable errors caused by the sky changes be-
tween measurements and imperfect alignment of rotating
optical elements.
The data presented in this paper were taken with the
R
-band filter. Magnitudes were calculated using calibrated
field stars either found in the literature or presented in PTF
(Palomar Transient Factory)
R
-band catalogue (
Ofek et al.
2012
) or USNO-B1.0 catalogue (
Monet et al. 2003
), depend-
ing on availability.
The exposure length was adjusted by the brightness of
each target, which was estimated during the short pointing
exposures, depending also on the sky conditions. The aver-
age photometric error in magnitudes is 0
.
04 mag. The data
were processed using the specialized pipeline described in
1
http://skinakas.physics.uoc.gr
detail by
King et al.
(
2014
) along with the telescope control
system.
Since we have introduced a Galactic latitude cut select-
ing objects with
|
b
|
>
10
◦
, the average colour excess in the
directions of our targets is relatively low,
E
(
B
−
V
) = 0
.
11
m
(
Schlafly & Finkbeiner 2011
), implying that the interstellar
polarization is less than 1
.
0% on average (
Serkowski et al.
1975
). The statistical uncertainty in the degree of polariza-
tion is less than 1% in most cases, while the EVPA is typi-
cally determined with a precision of 1
◦
– 10
◦
depending on
the source brightness and fractional polarization. Detail
ed
description of the instrument model and error analysis is
given in
King et al.
(
2014
).
In order to resolve the 180
◦
ambiguity of the EVPA we
followed a standard procedure (see e.g.
Abdo et al. 2010a
;
Ikejiri et al. 2011
;
Kiehlmann et al. 2013
), which is based
on the assumption that temporal variations of the EVPA
are smooth and gradual, hence adopting minimal changes of
the EVPA between consecutive measurements. We define the
EVPA variation as ∆
θ
n
=
|
θ
n+1
−
θ
n
|−
p
σ
(
θ
n+1
)
2
+
σ
(
θ
n
)
2
,
where
θ
n+1
and
θ
n
are the
n
+1 and
n
-th points of the EVPA
curve and
σ
(
θ
n+1
) and
σ
(
θ
n
) are the corresponding errors of
the position angles. If ∆
θ
n
>
90
◦
, we shift the angle
θ
n+1
by
±
n
×
180
◦
, where the integer
±
n
is chosen in such a
way that it minimizes ∆
θ
n
. If ∆
θ
n
≤
90
◦
, we leave
θ
n+1
unchanged.
Our first period of regular photometric and polarimetric
monitoring of blazars started in 2013 July and lasted until
the end of 2013 November. During the five-month period
we obtained more than 1100 measurements of 101 objects
from our sample almost uniformly spread over the observing
season of each object. The median cadence and total season
length for objects with ∆
t
-median smaller than 20 d (includ-
ing the June survey data,
Pavlidou et al. 2014
) is presented
in Fig.
1
, which is discussed in more detail in Sec.
3.3
.
2.3 Gamma-ray observations
The gamma-ray data were obtained with the Large Area
Telescope (LAT) onboard the
Fermi
gamma-ray space obser-
vatory, which observes the entire sky every 3 h at energies of
20 MeV – 300 GeV (
Atwood et al. 2009
). We analysed LAT
data in the energy range 100 MeV
≤
E
≤
100 GeV using the
unbinned likelihood analysis of the standard
Fermi
analy-
sis software package Science Tools v9r33p0 and the instru-
ment response function
P
7
REP
SOURCE
V
15. Source
class photons (evclass=2) were selected within a 15
◦
re-
gion of interest centred on a blazar. Cuts on the satel-
lite zenith angle (
<
100
◦
) and rocking angle (
<
52
◦
) were
used to exclude the Earth limb background. The diffuse
emission from the Galaxy was modelled using the spatial
model
gll
iem
v
05
rev
1. The extragalactic diffuse and resid-
ual instrumental backgrounds were included in the fit as an
isotropic spectral template
iso
source
v
05. The background
models
2
include all sources from the 2FGL catalogue within
15
◦
of the blazar. Photon fluxes of sources beyond 10
◦
from
the blazar and spectral shapes of all targets were fixed to
their values reported in 2FGL. The source is considered to
2
http://fermi.gsfc.nasa.gov/ssc/data/access/lat/
2yr_catalog/gll_psc_v07.xml
MNRAS
000
,
1
–
16
(2015)
4
D. Blinov et al.
Table 2.
Observational data for EVPA rotations detected by RoboPol i
n 2013. Columns (1),(2) - blazar identifiers; (3) - redshift;
(4)
- observational season length; (5) - average rotation rate;
(6) - total amplitude of EVPA change; (7) - number of observat
ions during
rotation; (8) - time duration of the rotation; (9) - TeV emiss
ion flag according to TeVCat
14
(“Y” means that the blazar has been detected
in gamma rays with
E >
1 TeV, “N” - otherwise); (10) - blazar subclass (LBL, IBL, HBL
denote low, intermediate and high synchrotron
peaked BL Lacertae objects, FSRQ – flat spectrum radio quasar
).
Blazar ID
Survey
z
T
obs
h
∆
θ
∆
t
i
∆
θ
max
N
points
T
rot
TeV
Class
name
(d) (deg/d) (deg)
(d)
RBPLJ0136+4751
OC 457
0.859
1
59
-6.6
-225
6
34
N FSRQ
13
RBPLJ0259+0747
PKS0256+075 0.893
2
72
-4.8
-180
6
38
N FSRQ
13
RBPLJ0721+7120*
S5 0716+71
0.31
3
88
-14.8
-208
11
14
Y
LBL
10
RBPLJ0854+2006*
OJ 287
0.306
4
51
-6.7
-154
10
23
N
LBL
10
RBPLJ1048+7143
S5 1044+71
1.15
5
142
-9.0
-188
22
21
N
−
RBPLJ1555+1111
PG 1553+113
−
129
5.6
128
8
23
Y
HBL
10
RBPLJ1558+5625
TXS1557+565
0.3
6
137
7.2
222
9
31
N
IBL?
11
RBPLJ1806+6949
3C 371
0.05
7
143
-16.5
-347
7
21
N
LBL
11
RBPLJ1806+6949
−′′−
−′′−
−′′−
13.3
238
5
18
N
−′′−
RBPLJ1927+6117
S4 1926+61
−
135
-4.4
-105
6
24
N
LBL
13
RBPLJ2202+4216
BL Lac
0.069
8
137
-51.0
-253
5
5
Y
LBL
10
RBPLJ2232+1143
CTA 102
1.037
1
140
-15.6
-312
8
20
N FSRQ
13
RBPLJ2232+1143
−′′−
−′′−
−′′−
-11.8
-140
6
12
N
−′′−
RBPLJ2243+2021 RGB J2243+203
−
169
-5.9
-183
5
31
N
LBL
12
RBPLJ2253+1608
3C 454.3
0.859
1
159
-18.3
-129
4
7
N FSRQ
13
RBPLJ2311+3425
B2 2308+34
1.817
9
36
3.3
74
20
23
N FSRQ
13
* Source belongs to sample (iii)
1
(
Hewitt & Burbidge 1987
);
2
(
Murphy et al. 1993
);
3
(
Nilsson et al. 2008
);
4
(
Nilsson et al. 2010
);
5
(
Polatidis et al. 1995
);
6
(
Falco et al. 1998
);
7
(
de Grijp et al. 1992
);
8
(
Vermeulen et al. 1995
);
9
(
Wills & Wills 1976
);
10
(
Donato et al. 2001
);
11
(
Ghisellini et al. 2011
);
12
(
Nieppola et al. 2006
);
13
(
Fan et al. 2012
);
14
http://tevcat.uchicago.edu
be detected if the test statistic, TS, provided by the analys
is
exceeds 10, which corresponds to approximately a 3
σ
detec-
tion level (
Nolan et al. 2012
). The systematic uncertainties
in the effective LAT area do not exceed 10 per cent in the
energy range we use (
Ackermann et al. 2012
). This makes
them insignificant with respect to the statistical errors, t
hat
dominate over the short time-scales analysed in this paper.
Moreover our analysis is based on the relative flux varia-
tions. Therefore the systematic uncertainties were not tak
en
into account.
Different time bins
t
int
, from 1 week to 25 d were used,
depending on the flux density of the object. In order to make
the analysis more robust we increased sampling of the pho-
ton flux curves shifting centres of the time bins by
t
int
/
4
interval from each other. This prevents losses of possible
short-term events in the light curves and reduces the depen-
dence of results on the particular position of the time bins.
The oversampling introduces an autocorrelation in the pho-
ton flux curves, which is however inessential for the analysi
s
used in this work.
3 RESULTS
3.1 Detected rotations of EVPA
The optical emission polarization plane of blazars is often
variable even within the course of a single night. There is no
objective physical definition of an EVPA rotation. Strictly
speaking, any change of the EVPA between two measure-
ments constitutes a rotation. However typically only high-
amplitude (
>
90
◦
), smooth and well tracked variations of
the EVPA are considered as rotations in the literature.
We accept a swing between two consecutive EVPA
measurements ∆
θ
=
|
θ
i+1
−
θ
i
|
as significant if ∆
θ >
p
σ
(
θ
i+1
)
2
+
σ
(
θ
i
)
2
. We define as an EVPA rotation any
continuous change of the EVPA curve with a total ampli-
tude ∆
θ
max
>
90
◦
, which is comprised by at least four mea-
surements with significant swings between them. Start and
end points of a rotation event are defined by a change of the
EVPA curve slope ∆
θ
i
/
∆
t
i
by a factor of 5 or a change of its
sign. This definition is rather conservative, and is in gener
al
consistent with rotations reported in the literature.
Using this definition, we identified 14 rotations of the
EVPA in 12 blazars from the main sample during the sea-
son of 2013 (see Table
2
). This number is comparable to the
number of previously known events of this type. Two more
blazars with detected rotations, namely RBPLJ0721+7120
and RBPLJ0854+2006, belong to the additional sample of
hand-picked sources. These blazars/events were not includ
ed
in the statistical or frequency analysis of the following se
c-
tions in this paper. The full season EVPA curves along
with the evolution of the polarization degree and the
R
-
band flux density, for all 14 blazars with detected rotations
,
are shown in Fig.
2
and listed in Table
2
. The EVPA rota-
tions are marked by filled black points. Clearly the events
we have considered as rotations based on our criteria are the
largest ∆
θ
max
rotation events that appear in these data sets.
They are all characterized by smooth variations with a well-
defined trend. Two events plotted in Fig.
2
do not follow the
definition strictly. These are the rotation events detected
in
the data sets of RBPLJ1048+7143 and RBPLJ2311+3425.
In both cases the rotations were interrupted by short, low
amplitude, albeit significant swings in the opposite direct
ion
with respect to the overall rotation. Since both events are
MNRAS
000
,
1
–
16
(2015)
RoboPol: EVPA rotations in blazars
5
Figure 2.
Evolution of polarization degree, polarization position a
ngle and
R
-band magnitude for blazars with a detected rotation in
the first RoboPol season. Periods of rotations are marked by fi
lled black points.
MNRAS
000
,
1
–
16
(2015)
6
D. Blinov et al.
Figure 2
–
continued
(Continued) Evolution of polarization degree, polarizati
on position angle and
R
-band magnitude for blazars with
a detected rotation in the first RoboPol season. Periods of ro
tations are marked by filled black points.
MNRAS
000
,
1
–
16
(2015)
RoboPol: EVPA rotations in blazars
7
Figure 2
–
continued
(Continued) Evolution of polarization degree, polarizati
on position angle and
R
-band magnitude for blazars with
a detected rotation in the first RoboPol season. Periods of ro
tations are marked by filled black points.
well sampled these small deviations do not introduce any
significant difference in the overall EVPA trend. Hence both
events can be considered as single, large ∆
θ
max
rotations. In
addition the RBPLJ2311+3425 event has an amplitude of
∼
74
◦
, which is less than the lower limit we accepted. How-
ever the start and end points of the rotation are not defined
due to a sparse sampling. It is likely that this well defined
EVPA change would meet the 90
◦
limit if we had a longer
data set for this object. It is for this reason that we include
this event in our sample of rotations. Both events have not
been used in any of our statistical analyses involving com-
parison between simulated and observed rotations.
Some of the EVPA rotation events are coincident with
an increase in the total flux, as it follows from a visual in-
spection of Fig.
2
. A quantitative comparison between the
optical flux and the polarization variations will be present
ed
in a forthcoming paper.
3.2 General properties of EVPA rotations and
rotators
We estimated the maximal amplitude ∆
θ
max
and the dura-
tion of the rotations
T
rot
, using the first and last points of
each event. Due to a moderate sampling and 180
◦
EVPA am-
biguity, the rotation start and/or end points cannot be pin-
pointed accurately in five cases (namely RBPLJ0136+4761,
RBPLJ0259+0747, RBPLJ1048+7143, RBPLJ1806+6949
and RBPLJ2311+3425). This ambiguity affects the esti-
mated ∆
θ
max
and
T
rot
of the event, which should really be
considered as lower limits in this case. We also estimated
the average rotation rate as
h
∆
θ/
∆
t
i
= ∆
θ
max
/T
rot
. These
parameters as well as the blazar class and the TeV emission
flag are listed in Table
2
.
We also collected data from the literature on previously
known rotations of EVPA in blazars which show this be-
haviour (“rotators” hereafter). Rates and ∆
θ
max
of these ro-
tations were estimated from plots in the respective papers.
These parameters as well as the blazar class and the TeV
emission flag are listed in Table
3
.
Figure 3.
Distributions of amplitudes and rates of EVPA ro-
tations detected in RoboPol’s first season and reported in th
e
literature.
The distribution of ∆
θ
max
and rates of EVPA rotations
from historical and RoboPol data are shown in Fig.
3
. The
number of detected rotations clearly decreases with grow-
ing ∆
θ
max
. At the same time slow rotations dominate in the
sample. This is presumably caused by a selection effect, be-
cause fast rotations require better sampling of observatio
ns.
Summarizing data on all known EVPA rotations in
blazars to date we can list the following properties:
(i) all known blazars with detected EVPA rotations are
in the 2FGL catalogue (i.e. they are “gamma-ray–loud”
sources);
(ii) there are blazars known as TeV emitters as well as
non-TeV sources among rotators;
(iii) all subclasses of blazars show rotations of the EVPA,
regardless of the position of the synchrotron peak maximum
or the BL Lac/FSRQ dichotomy;
(iv) there are eight blazars with more than one rotation
detected. Comparison of these rotations shows that a single
source can show rotations in both directions (five blazars
known so far with this behaviour) and rotations observed
in the same source can be of significantly different rates
MNRAS
000
,
1
–
16
(2015)
8
D. Blinov et al.
Table 3.
Data on known rotation of optical EVPA in blazars. Columns (1
),(2) - blazar identifiers; (3) -
average rotation rate; (4) - total amplitude of EVPA change;
(5) - TeV emission flag according to TeVCat
5
(“Y” means that the blazar has been detected in gamma rays wit
h
E >
1 TeV, “N” - otherwise); (6) - blazar
subclass (LBL, IBL, HBL denote low, intermediate and high sy
nchrotron peaked BL Lacertae objects, FSRQ
– flat spectrum radio quasar); (7) - reference.
Blazar ID
Survey
h
∆
θ
∆
t
i
∆
θ
max
TeV
Class
Reference
name
(deg d
−
1
) (deg)
RBPLJ0423
−
0120 PKS 0420
−
014
−
11.1
−
110
N
FSRQ
4
(
D’Arcangelo et al. 2007
)
RBPLJ0721+7120
S5 0716+71
+130
+180
Y
LBL
1
(
Larionov et al. 2013b
)
RBPLJ0854+2006
OJ 287
−
17
−
180
N
LBL
1
(
Kikuchi et al. 1988
)
RBPLJ0958+6533
S4 0954+65
+18.2
+240
N
LBL
2
(
Larionov et al. 2011
)
RBPLJ1221+2813
W Comae
≥
+3
.
0
+110
Y
IBL
3
(
Ben ́ıtez et al. 2013
)
RBPLJ1256
−
0547
3C 279
−
9
−
180
Y
FSRQ
4
(
Abdo et al. 2010a
)
RBPLJ1256
−
0547
3C 279
+4.3
+290
−′′−
−′′−
(
Larionov et al. 2008
)
RBPLJ1256
−
0547
3C 279
+4.7
+140
−′′−
−′′−
(
Aleksi ́c et al. 2014a
)
RBPLJ1512
−
0905 PKS 1510
−
089
+15.6
+720
Y
FSRQ
4
(
Marscher et al. 2010
)
RBPLJ1512
−
0905 PKS 1510
−
089
+12
+400
−′′−
−′′−
(
Aleksi ́c et al. 2014b
)
RBPLJ1512
−
0905 PKS 1510
−
089
−
50
−
250
−′′−
−′′−
(
Aleksi ́c et al. 2014b
)
RBPLJ1512
−
0905 PKS 1510
−
089
+11.7
+500
−′′−
−′′−
(
Sasada et al. 2011
)
6
RBPLJ2202+4216
BL Lac
+46
+220
Y
IBL
1
(
Marscher et al. 2008
)
RBPLJ2202+4216
BL Lac
+21
+210
−′′−
−′′−
(
Sillanp
̈
a
̈
a et al. 1993
)
RBPLJ2232+1143
CTA 102
−
60
−
180
N
FSRQ
4
(
Larionov et al. 2013a
)
RBPLJ2253+1608
3C 454.3
+16.3
+130
N
FSRQ
4
(
Sasada et al. 2010
)
RBPLJ2253+1608
3C 454.3
+9.3
+400
−′′−
−′′−
(
Sasada et al. 2012
)
1
(
Donato et al. 2001
);
2
(
Ghisellini et al. 2011
);
3
(
Tagliaferri et al. 2000
);
4
(
Fan et al. 2012
)
5
http://tevcat.uchicago.edu
;
6
same as in
Marscher et al.
(
2010
).
(in seven blazars rates differ by a factor larger than two
in speed).
3.3 Observed frequency of EVPA rotations
The efficiency of an EVPA rotation detection depends on
the intrinsic rate of the rotation as well as the frequency
and uniformity of the observing cadence. The ambiguity of
the polarization position angle introduces an upper limit o
n
the rotation rate that can be unequivocally detected with a
given typical cadence of observations. Clearly, for a typic
al
time interval between observations
h
∆
t
i
, no EVPA rotation
with a rate higher than 90
◦
/
h
∆
t
i
can be observed.
For each blazar in our sample we found the median time
difference between successive observations ∆
t
-median and
the total observing season length (defined as the time differ-
ence between the first and the last observations)
T
obs
. These
quantities (for blazars observed with ∆
t
-median
≤
20 days)
are shown in Fig.
1
. In the same figure, we also plot three
lines which indicate the necessary ∆
t
-median and
T
obs
for
detection of EVPA rotations at rates
≤
10 (solid line),
≤
15
(dashed line) and
≤
20 (dotted line) degrees per day.
The leftmost vertical part of each line represents the
shortest
T
obs
needed to detect a rotation of ∆
θ
max
= 90
◦
at a given rotation rate. The inclined portion of each line
is determined by our requirement on a rotation event to be
comprised by a minimum of four points. Given this require-
ment, as ∆
t
-median increases, so does
T
obs
. An EVPA data
set with ∆
t
-median and
T
obs
on that line can allow detec-
tion of EVPA rotations with 90
◦
≤
∆
θ
max
≤
270
◦
. The
horizontal part indicates the maximum ∆
t
-median allowed
the detection of a rotation event under the requirement of
∆
θ
≤
90
◦
in EVPA between two consecutive points.
We can now estimate the frequency with which EVPA
rotations appear in blazars as follows. Out of the 14 detecte
d
rotations in blazars of the main sample, 8 have rates less
than 10 deg d
−
1
. There are also 41 main sample (“gamma-
ray bright”) blazars that were observed with ∆
t
-median and
T
obs
(see Table
4
) within the region defined by the solid
line in Fig.
1
. The total observing length for these blazars is
6432 d. Thereby we estimate the frequency of “slow” rota-
tions (rate
<
10 deg d
−
1
) in the main sample sources as one
rotation in
∼
800 days (6432 d / 8 rotations). Following the
same reasoning we estimate average frequencies of rotation
s
for blazars in the main sample with rates
<
15 deg d
−
1
and
<
20 deg d
−
1
as one rotation in
∼
490 d (4912/10) and
∼
180 d (2363/13), respectively.
3.4 EVPA variability in blazars of different
samples
In order to address the question whether “the EVPA vari-
ability is different in objects where rotations were detecte
d
compared to the rest of the main sample and to the control
sample” we collate all EVPA “swing” events and measure
their ∆
θ
max
and rates. We define an EVPA “swing” as any
continuous change of the EVPA curve, without a lower limit
in its ∆
θ
max
or in the number of measurements. As before,
start and end points of a swing event are defined by a change
of the EVPA curve slope by a factor of 5 or a change of its
sign.
We identified all such events for all blazars of the main
and control samples within the 10 deg d
−
1
“detection box”in
Fig.
1
, and measured their amplitude, ∆
θ
max
, and mean ro-
tation rate. The cumulative distribution function (hereaf
ter
CDF; e.g.
Wall & Jenkins 2012
) of the EVPA swings ∆
θ
max
and rotation rates for blazars in the main sample which
showed rotations (“rotators”), blazars in the main sample,
MNRAS
000
,
1
–
16
(2015)
RoboPol: EVPA rotations in blazars
9
Table 4.
Sources of the main and control samples within the
h
∆
θ
∆
t
i
<
10 deg d
−
1
“detection box”.
Blazar ID
Survey
T
obs
h
∆
t
i
Blazar ID
Survey
T
obs
h
∆
t
i
(RBPL. . . )
name
(d)
(d) (RBPL. . . )
name
(d)
(d)
Main sample
J1838+4802
GB6 J1838+4802
121 7.0
J0045+2127
GB6 J0045+2127
33
4.0 J1841+3218
RXJ1841.7+321
8
152 6.0
J0114+1325
GB6 J0114+1325
38
6.0 J1903+5540
TXS 1902+556
13
5 5.0
J0211+1051 MG1 J021114+1051 85
4.0 J1959+6508
1ES 1959+650
143 5.0
J0217+0837
ZS0214+083
85
6.0 J2005+7752
S5 2007+77
140 7.0
J0423
−
0120
PKS 0420
−
01
12
4.0 J2015
−
0137
PKS2012
−
017
155 6.5
J0841+7053
4C 71.07
71
6.0 J2016
−
0903
PMNJ2016
−
0903
155 7.0
J1512
−
0905
PKS 1510
−
08
88
2.0 J2022+7611
S5 2023+760
158 7.0
J1542+6129
GB6 J1542+6129
87
4.0 J2030
−
0622
TXS 2027
−
065
143 5.0
J1553+1256
PKS 1551+130
132 4.0 J2039
−
1046
TXS 2036
−
109
144 5.5
J1604+5714
GB6 J1604+5714
135 7.0 J2131
−
0915
RBS1752
127 5.0
J1607+1551
4C 15.54
136 8.0 J2143+1743
OX 169
119 5.0
J1635+3808
4C 38.41
121 2.0 J2148+0657
4C 6.69
152 4.5
J1642+3948
3C 345
148 6.0 J2149+0322
PKSB 2147+031
169 6.5
J1653+3945
Mkn 501
153 4.0 J2150
−
1410
TXS 2147
−
144
130 8.0
J1725+1152
1H 1720+117
120 3.0 J2225
−
0457
3C 446
144 6.0
J1748+7005
S4 1749+70
87
3.0 J2251+4030 MG4 J225201+4030 17
7 6.5
J1751+0939
OT 081
154 4.0 J2334+0736
TXS 2331+073
138 8.0
J1754+3212
RXJ1754.1+3212
134 5.0 J2340+8015
BZBJ2340+80
15
113 5.5
J1800+7828
S5 1803+784
133 5.0
Control sample
J1809+2041
RXJ1809.3+2041
152 4.5 J1551+5806
SBS1550+582
118 5.0
J1813+3144
B2 1811+31
150 6.0 J1638+5720
S4 1637+57
138 4.5
J1836+3136
RXJ1836.2+3136
151 6.0 J2042+7508
4C +74.26
99
5
.0
which did not show rotations (“non-rotators”), as well as fo
r
blazars in the control sample, are shown in Fig.
4
.
We performed a two sample Kolmogorov–Smirnov (K–
S) test (e.g.
Wall & Jenkins 2012
) pairwise for three sam-
ples of collected swing amplitudes and rates with the null-
hypothesis that these samples are drawn from the same dis-
tribution. The null-hypothesis is rejected for rotators an
d
non-rotators with the
p
-value = 1
.
2
×
10
−
5
, and for rota-
tors and the control sample (
p
-value = 5
×
10
−
3
). At the
same time the distribution of swing amplitudes in the non-
rotators and control sample sources is indistinguishable a
c-
cording to the test (
p
-value = 0
.
35). The maximum differ-
ence between the CDFs of non-rotators and rotators is 0.29.
It is reached at ∆
θ
max
≈
25
◦
. Even if we exclude the 14
rotations (i.e. the largest ∆
θ
max
swings) of the main sample
blazars, rotators still remain different from the non-rotat
ors
(
p
-value = 2
×
10
−
3
).
A similar analysis (as the one for ∆
θ
max
) for the distri-
butions of EVPA swing rates leads to the same conclusion.
The null-hypothesis is rejected for the rotators and the non
-
rotators (
p
-value = 1
.
4
×
10
−
6
) and rotators vs. the control
sample (
p
-value = 5
×
10
−
3
), while it can not be rejected for
the non-rotators and control sample (
p
-value = 0
.
18).
We therefore conclude that blazars with detected rota-
tions show significantly larger
∆
θ
max
and faster EVPA vari-
ations when compared to blazars with no detected rotations.
This difference cannot be attributed to differences in the
sampling properties of the data sets. Therefore, the lack of
detection of EVPA rotations in the “non-rotators” member
of the main sample, as well as the blazar in the control sam-
ple, may have a physical origin. Most of the non-rotators
in the main and control samples may never show an EVPA
rotation.
4 RANDOM WALKS AS THE ORIGIN OF
EVPA ROTATIONS
4.1 MC simulations of EVPA swings
Potentially EVPA swings can be explained by a stochastic
process, which is physically justified by a presence of many
independent cells in the emission region (e.g.
Jones et al.
1985
;
D’Arcangelo et al. 2007
). According to this interpre-
tation, the magnetic field is turbulent and apparent rotatio
ns
result from a random walk of the full polarization vector di-
rection as new cells with random magnetic field orientations
appear in the emission region (
Marscher 2014
). In order to
estimate the probabilities that the EVPA rotations we ob-
served with RoboPol are produced by this kind of multicell
random walk process we performed MC simulations of the
stochastic variability of the polarization vector on the QU
plane following
Kiehlmann et al.
(
2013
).
For each blazar where an EVPA rotation event was
observed, we created 10
4
artificial light curves, each one
with duration
T
obs
. The time steps ∆
t
i
between consecutive
points were drawn from a truncated power-law distribution,
which approximates well the distribution of the time steps i
n
all observed lightcurves. The parameters of this distribut
ion
(∆
t
min
, ∆
t
max
and the power-law index) were determined by
fitting it to the distribution of observed ∆
t
i
for each object.
The total flux density
I
i
emitted at each time step ∆
t
i
,
was drawn from a log-normal distribution. Such a distribu-
tion approximates reasonably well the distribution of the o
b-
served flux densities for all blazars. The mean and variance
of the log-normal distribution was set equal to the sample
mean and variance of the distribution of the flux density of
each blazar.
The maximum possible fractional polarization produced
by a uniform magnetic field is
P
max
= (
α
+ 1)
/
(
α
+ 5
/
3)
≈
MNRAS
000
,
1
–
16
(2015)
10
D. Blinov et al.
Figure 4.
CDFs of ∆
θ
max
and average rates for the main sample
rotators (95 EVPA swings), the main sample non-rotators (29
8
swings) and control sample sources (11 swings). See Sec.
3.4
for
details.
0
.
78 (
Pacholczyk 1970
). In the case of unresolved emis-
sion region comprising
N
independent cells with a uni-
form magnetic field, but randomly oriented among them,
the average fractional polarization is given by the equatio
n
(
Hughes & Miller 1991
):
h
P
obs
i ≈
P
max
√
N
.
(1)
We used this equation and the observed average polarization
fraction,
h
P
obs
i
, to estimate the number of cells,
N
, for each
blazar. Each
k
-th cell at
i
-th time step was assigned a flux
density
I
i
,
k
(which was set equal to
I
i
/N
for all cells at each
time step) and a set of fractional Stokes parameters
q
i
,
k
and
u
i
,
k
. They were found as
q
i
,
k
=
q
0
i
,
k
P
max
√
(
q
0
i
,
k
)
2
+(
u
0
i
,
k
)
2
u
i
,
k
=
u
0
i
,
k
P
max
√
(
q
0
i
,
k
)
2
+(
u
0
i
,
k
)
2
,
(2)
where
q
0
i
,
k
and
u
0
i
,
k
are two numbers drawn from the standard
normal distribution. Thereby the emission of each cell has
polarization fraction
P
max
. The sums
Q
i
=
I
i
P
N
k=1
q
i
,
k
and
U
i
=
I
i
P
N
k=1
u
i
,
k
determine the total Stokes parameters of
the emitting region at each time step.
At each time step the Stokes parameters of
N
var
(∆
t
i
)
cells, selected randomly, were replaced by new values. The
number of cells for replacement was estimated (from the
average variance of the polarization degree) as follows:
N
var
(∆
t
i
) =
∆
t
i
h
∆
t
i
σ
(
P
obs
)
h
P
obs
i
N,
(3)
where
σ
(
P
obs
) is the observed standard deviation of the de-
gree of polarization for each blazar, and
h
∆
t
i
is the average
time difference between observations.
It was confirmed that the simulated and observational
data in corresponding blazars have similar statistical pro
p-
erties. Namely, the standard deviation and average of the
polarization fraction are consistent with
σ
(
P
obs
) and
h
P
obs
i
.
4.1.1 Individual rotations
Using the algorithm described in Sec.
3
, i.e. the same algo-
rithm we used to identify rotations in real data, we iden-
tified all rotations in the simulated data and found the
number
N
rot
of “successful” data sets, where at least one
Table 5.
Random walk modelling results for EVPA rotations
detected by RoboPol in 2013. (1) - blazar identifier; (2) - occ
ur-
rence of rotations with ∆
θ
max
,
simul
≥
∆
θ
max
,
obs
estimated from
the simulations; (3) - probability that a rotation produced
by the
random walk will be observed in
T
obs
.
Blazar ID
T
occ
P(RW)
(d)
RBPLJ0136+4751 505
0.11
RBPLJ0259+0747 151
0.48
RBPLJ0721+7120 325
0.28
RBPLJ0854+2006 142
0.36
RBPLJ1048+7143 180
0.79
RBPLJ1555+1111 128
1.00
RBPLJ1558+5625 266
0.51
RBPLJ1806+6949 965
0.15
RBPLJ1806+6949 259
0.55
RBPLJ1927+6117 137
0.98
RBPLJ2202+4216 633
0.21
RBPLJ2232+1143 1557
0.09
RBPLJ2232+1143 178
0.87
RBPLJ2243+2021 183
0.92
RBPLJ2253+1608 184
0.86
RBPLJ2311+3425 61
0.74
rotation with ∆
θ
max
larger or equal to ∆
θ
max
,
obs
was de-
tected. We then estimated two ratios:
P
(RW) = N
rot
/
10
4
and
T
occ
= 10
4
·
T
obs
/N
rot
. The first ratio determines the
probability to observe an EVPA rotation due to a random
walk for each one of the observed EVPA curves for the given
T
obs
. The second ratio determines the average time interval
between random walk rotations (i.e. the average occurrence
rate for each blazar). The probabilities
P
(RW) and
T
occ
are
listed in Table
5
. The probabilities are larger than 10% in
all but one object, and in some cases, they approach unity.
This result indicates that the rotations we observed in some
objects could be the result of a random walk process.
4.1.2 Rotations as a population
In this section we test the hypothesis that
all
the rotations
observed by RoboPol in blazars of the main sample are pro-
duced by the stochastic process. According to the analysis
in Sec.
3.4
blazars exhibiting rotations have different prop-
erties when compared to non-rotators. Therefore the sample
of rotators must be considered separately.
We performed the following simulation. At each itera-
tion, an artificial EVPA curve was generated individually fo
r
each rotator from the main sample as explained in Sec.
4.1
.
In each of the simuated EVPA curves we identified the
largest rotation and constructed the CDF of ∆
θ
max
,
simul
among the blazars. An iteration was considered to be “suc-
cessful” only in the case when the CDF of ∆
θ
max
,
simul
was
lower or equal to the CDF of ∆
θ
max
,
obs
, i.e. the simulated set
of EVPA curves had higher or equal fraction of rotations of
a given length compared to the observed set. In the cases of
RBPLJ1806+6949 and RBPLJ2232+1143 where double ro-
tations were observed, we simulated only the largest ∆
θ
max
rotations.
The CDF of ∆
θ
max
,
obs
along with a subset of 100 simu-
lated CDFs is shown in Fig.
5
. It was found that only 1
.
5% in
10
4
trials were “successful”. Therefore, the probability that
MNRAS
000
,
1
–
16
(2015)