ar X iv :1 80 6. 05 36 7v 2 [ as tr o- ph .G A ] 3 S ep 2 01 8 Astronomy & Astrophysics manuscript no. pks1510-low c©ESO 2018 September 5, 2018 Detection of persistent VHE gamma-ray emission from PKS 1510–089 by the MAGIC telescopes during low states between 2012 and 2017. MAGIC Collaboration: V. A. Acciari1, S. Ansoldi2, 21, L. A. Antonelli3, A. Arbet Engels4, C. Arcaro5, D. Baack6, A. Babić7, B. Banerjee8, P. Bangale9, U. Barres de Almeida9, 10, J. A. Barrio11, W. Bednarek12, E. Bernardini5, 13, 23, A. Berti2,24, J. Besenrieder9, W. Bhattacharyya13 , C. Bigongiari3, A. Biland4, O. Blanch14, G. Bonnoli15, R. Carosi16, G. Ceribella9, S. Cikota7, S. M. Colak14, P. Colin9, E. Colombo1, J. L. Contreras11, J. Cortina14, S. Covino3, V. D’Elia3, P. Da Vela16, F. Dazzi3, A. De Angelis5, B. De Lotto2, M. Delfino14, 25, J. Delgado14, 25, F. Di Pierro, E. Do Souto Espiñera14, A. Domínguez11, D. Dominis Prester7, D. Dorner17, M. Doro5, S. Einecke6, D. Elsaesser6, V. Fallah Ramazani18, A. Fattorini6, A. Fernández-Barral5, 14, G. Ferrara3, D. Fidalgo11, L. Foffano5, M. V. Fonseca11, L. Font19, C. Fruck9, D. Galindo20, S. Gallozzi3, R. J. García López1, M. Garczarczyk13, M. Gaug19, P. Giammaria3, N. Godinović7, D. Guberman14, D. Hadasch21, A. Hahn9, T. Hassan14, J. Herrera1, J. Hoang11, D. Hrupec7, S. Inoue21, K. Ishio9, Y. Iwamura21, H. Kubo21, J. Kushida21, D. Kuveždić7, A. Lamastra3, D. Lelas7, F. Leone3, E. Lindfors18, S. Lombardi3, F. Longo2,24, M. López11, A. López-Oramas1, C. Maggio19, P. Majumdar8, M. Makariev22, G. Maneva22, M. Manganaro7, K. Mannheim17, L. Maraschi3, M. Mariotti5, M. Martínez14, S. Masuda21, D. Mazin9, 21, M. Minev22, J. M. Miranda15, R. Mirzoyan9, E. Molina20, A. Moralejo14, V. Moreno19, E. Moretti14, P. Munar-Adrover19, V. Neustroev18, A. Niedzwiecki12, M. Nievas Rosillo11, C. Nigro13⋆, K. Nilsson18, D. Ninci14, K. Nishijima21, K. Noda21, L. Nogués14, S. Paiano5, J. Palacio14, D. Paneque9, R. Paoletti15, J. M. Paredes20, G. Pedaletti13, P. Peñil11, M. Peresano2, M. Persic2, 26, P. G. Prada Moroni16, E. Prandini5, I. Puljak7, J. R. Garcia9, W. Rhode6, M. Ribó20, J. Rico14, C. Righi3, A. Rugliancich16 , L. Saha11, T. Saito21, K. Satalecka13, T. Schweizer9, J. Sitarek12⋆, I. Šnidarić7, D. Sobczynska12 , A. Somero1, A. Stamerra3, M. Strzys9, T. Surić7, F. Tavecchio3, P. Temnikov22, T. Terzić7, M. Teshima9, 21, N. Torres-Albà20, S. Tsujimoto21, J. van Scherpenberg9 , G. Vanzo1, M. Vazquez Acosta1, I. Vovk9, J. E. Ward14, M. Will9, D. Zarić7; and Fermi-LAT Collaboration: J. Becerra González1 ⋆; and C. M. Raiteri27, A. Sandrinelli28, 29, T. Hovatta31, S. Kiehlmann30, W. Max-Moerbeck32 , M. Tornikoski33, A. Lähteenmäki33, 34, 35, J. Tammi33, V. Ramakrishnan33, C. Thum36, I. Agudo37, S. N. Molina37, J. L. Gómez37, A. Fuentes37, C. Casadio38, E. Traianou38, I. Myserlis38, and J.-Y. Kim38 (Affiliations can be found after the references) Received ; accepted ABSTRACT Context. PKS 1510–089 is a flat spectrum radio quasar strongly variable in the optical and GeV range. So far, very-high-energy (VHE, > 100 GeV) emission has been observed from this source during either long high states of optical and GeV activity or during short flares. Aims. We search for low-state VHE gamma-ray emission from PKS 1510–089. We aim to characterize and model the source in a broad-band context, which would provide a baseline over which high states and flares could be better understood. Methods. PKS 1510–089 has been monitored by the MAGIC telescopes since 2012. We use daily binned Fermi-LAT flux measurements of PKS 1510–089 to characterize the GeV emission and select the observation periods of MAGIC during low state of activity. For the selected times we compute the average radio, IR, optical, UV, X-ray and gamma-ray emission to construct a low-state spectral energy distribution of the source. The broadband emission is modelled within an External Compton scenario with a stationary emission region through which plasma and magnetic field are flowing. We perform also the emission-model-independent calculations of the maximum absorption in the broad line region (BLR) using two different models. Results. The MAGIC telescopes collected 75 hrs of data during times when the Fermi-LAT flux measured above 1 GeV was below 3×10−8cm−2s−1, which is the threshold adopted for the definition of a low gamma-ray activity state. The data show a strongly significant (9.5σ) VHE gamma-ray emission at the level of (4.27 ± 0.61stat) × 10−12 cm−2s−1 above 150 GeV, a factor 80 smaller than the highest flare observed so far from this object. Despite the lower flux, the spectral shape is consistent with earlier detections in the VHE band. The broad-band emission is compatible with the External Compton scenario assuming a large emission region located beyond the broad line region. For the first time the gamma-ray data allow us to place a limit on the location of the emission region during a low gamma-ray state of a FSRQ. For the used model of the BLR, the 95% C.L. on the location of the emission region allows us to place it at the distance > 74% of the outer radius of the BLR. Key words. galaxies: active – galaxies: jets – gamma rays: galaxies – quasars: individual: PKS 1510-089 Article number, page 1 of 13 http://arxiv.org/abs/1806.05367v2 A&A proofs: manuscript no. pks1510-low Article number, page 2 of 13 MAGIC Collaboration: V. A. Acciari et al.: Detection of VHE gamma rays from PKS 1510–089 during low state. 1. Introduction PKS 1510–089 is a bright flat spectrum radio quasar (FSRQ) located at a redshift of z = 0.36 (Tanner et al. 1996). It was the second FSRQ to be detected in the very-high- energy (VHE, > 100 GeV) range (Abramowski et al. 2013). The source is monitored by various instruments spanning the full range from radio up to VHE gamma rays (see e.g. Marscher et al. 2010; Aleksić et al. 2014; Ahnen et al. 2017a). Similarly to other FSRQs, the GeV gamma-ray emission of PKS 1510–089 is strongly variable (Abdo et al. 2010; Brown 2013; Saito et al. 2013; Prince et al. 2017). Multiple optical flares have been observed from PKS 1510–089 (Liller & Liller 1975; Zacharias et al. 2016)1. Significant VHE gamma-ray emission from PKS 1510– 089 has been observed on a few occasions: during en- hanced optical and GeV states in 2009 (Abramowski et al. 2013) and 2012 (Aleksić et al. 2014) and during short flares in 2015 (Ahnen et al. 2017a; Zacharias et al. 2016) and 2016 (Zacharias et al. 2017). Interestingly, no variability in VHE gamma rays has been observed during (or between) the high optical/GeV states in 2009 and 2012 (Abramowski et al. 2013; Aleksić et al. 2014). The GeV state of PKS 1510–089 can be studied using the Fermi-LAT all-sky monitoring data. MAGIC is a system of two Imaging Atmospheric Cherenkov Telescopes designed for ob- servations of gamma rays with energies from a few tens of GeV up to a few tens of TeV (Aleksić et al. 2016a). Since the detec- tion of VHE gamma-ray emission from PKS 1510–089 in 2012, a monitoring program is being performed with the MAGIC tele- scopes. The aimed cadence of monitoring is 2-6 pointings per month, with individual exposures of 1-3 hrs. The source is visi- ble by MAGIC 5 months of the year. We use the Fermi-LAT data to select periods of low gamma-ray emission of PKS 1510–089. Then, we select a subsample of the MAGIC telescope data taken between 2012 and 2017, and contemporaneous multiwavelength data from a number of other instruments, in order to study the quiescent VHE gamma-ray state of the source. Such low emis- sion can then be used as a baseline for modeling of flaring states. In Section 2 we briefly introduce the instruments that pro- vided multiwavelength data, describe the data reduction proce- dures and explain the principle of low-state data selection. In Section 3 we present the results of the observations, and the broadband emission modelling is illustrated in Section 4. The most important findings are summarized in Section 5. 2. Data The continuous monitoring of PKS 1510–089 in the GeV band provided by Fermi-LAT allows us to identify the low emission states of the source. Multiwavelength light curves from the radio band up to the GeV band are shown in Fig. 1. 2.1. Fermi-LAT Fermi-LAT monitors the high energy gamma-ray sky in the en- ergy range from 20 MeV to beyond 300 GeV (Atwood et al. 2009). For this work, we have analyzed the Pass 8 SOURCE class events within a region of interest (ROI) of 10◦ ra- ⋆ Corresponding authors: J. Sitarek (jsitarek@uni.lodz.pl), C. Nigro (cosimo.nigro@desy.de), J. Becerra Gonzalez (jbecerra@iac.es) 1 see also http://users.utu.fi/kani/1m/PKS_1510-089_jy. html dius centered at the position of PKS 1510–089 in the en- ergy range from 100 MeV to 300 GeV. A zenith angle cut of < 90◦ was applied to reduce the contamination from the Earth’s limb. The analysis was performed with the ScienceTools software package version v11r7p0 using the P8R2_SOURCE_V62 instrument response function and the gll_iem_v06 and iso_P8R2_SOURCE_V6_v06models3 for the Galactic and isotropic diffuse emission (Acero et al. 2016), re- spectively. A first unbinned likelihood fit was performed for the events collected within five months from 01 February to 30 June 2013 (MJD 56324–56474) using gtlike, and including in the model file all 3FGL sources (Acero et al. 2015) of 20◦ from PKS 1510– 089. We repeat the same 5-month analysis using the preliminary 8-year Source List 4 instead of the 3FGL catalog to search for bright sources within 20◦ of PKS 1510–089 . No new strong sources were found. The model generated from the 3FGL cat- alog was used for the subsequent analysis. As we are interested in short time scale (daily) fluxes of PKS 1510–089 the purpose of this first fit is to identify weak nearby sources that can be re- moved from the source model, simplifying it. Hence, the sources with a test statistic (TS; Mattox et al. 1996) below 5 were re- moved from the model file. Next, the optimized output model file was used to produce the PKS 1510–089 light curve with 1- day time bins above 1 GeV in the full time period from 5 De- cember 2011 to 7 August 2017 (MJD 55900–57972). The same optimized output model is later also used for the spectral analy- sis. In the light curve calculations the spectra of PKS 1510–089 were modeled as power law leaving both the flux normaliza- tion and the spectral index as free parameters. The normaliza- tion of the Galactic and isotropic diffuse emission models was left to vary freely during the calculation of both the light curves and the spectrum. In addition, the spectra of all sources except PKS 1510–089 and the highly variable source 3FGL 1532.7- 1319 (located 6.45◦ from PKS 1510–089, and having variability index of 1924.7 from the 3FGL catalog) were fixed to the catalog values. In order to estimate when the flux can be considered being in a low state, we first calculate a light curve in relatively wide bins of 30 days in the full time period. This allows us to estimate the flux with relative uncertainty . 20% for all the points and hence disentangle intrinsic variability from the fluctuations of the measured flux induced by statistical uncertainties. In Fig. 2 we present the distribution of the flux above 1 GeV, which shows that during the low state the flux is between (1–3)×10−8 cm−2s−1 in contrast to the value of > 3× 10−8cm−2s−1 during active (flar- ing) periods. Hence, we select the days of low state if: F>1GeV < 3 × 10−8cm−2s−1. (1) The cut separates the low flux peak of the daily fluxes distri- bution from the power-law extension of the high-flux days (see Fig. 2). The effect that choosing a different energy threshold would have on the data selection is discussed in Appendix A. We note, however, that due to the low state and short exposure times the flux measurements during single nights are quite un- certain. The typical uncertainty of the flux in those time bins 2 http://fermi.gsfc.nasa.gov/ssc/data/analysis/ documentation/Cicerone/Cicerone_LAT_IRFs/IRF_overview. html 3 http://fermi.gsfc.nasa.gov/ssc/data/access/lat/ BackgroundModels.html 4 https://fermi.gsfc.nasa.gov/ssc/data/access/lat/ fl8y/ Article number, page 3 of 13 http://users.utu.fi/kani/1m/PKS_1510-089_jy.html http://users.utu.fi/kani/1m/PKS_1510-089_jy.html http://fermi.gsfc.nasa.gov/ssc/data/analysis/documentation/Cicerone/Cicerone_LAT_IRFs/IRF_overview.html http://fermi.gsfc.nasa.gov/ssc/data/analysis/documentation/Cicerone/Cicerone_LAT_IRFs/IRF_overview.html http://fermi.gsfc.nasa.gov/ssc/data/analysis/documentation/Cicerone/Cicerone_LAT_IRFs/IRF_overview.html http://fermi.gsfc.nasa.gov/ssc/data/access/lat/BackgroundModels.html http://fermi.gsfc.nasa.gov/ssc/data/access/lat/BackgroundModels.html https://fermi.gsfc.nasa.gov/ssc/data/access/lat/fl8y/ https://fermi.gsfc.nasa.gov/ssc/data/access/lat/fl8y/ A&A proofs: manuscript no. pks1510-low 5580056000562005640056600568005700057200574005760057800 9−10 8−10 7−10 6−10 ] -1 s -2 F [c m Fermi >1 GeV 55800 56000 56200 56400 56600 56800 57000 57200 57400 57600 57800 0 5 10 15 20 ] -1 s -2 er g cm -1 2 F [1 0 Swift-XRT 2-10 keV 55800 56000 56200 56400 56600 56800 57000 57200 57400 57600 57800 0 2 4 6 F [m Jy ] UVOT U-band 55800 56000 56200 56400 56600 56800 57000 57200 57400 57600 57800 1 10 F [m Jy ] R-band (KVA, MAPCAT, SMARTS) 55800 56000 56200 56400 56600 56800 57000 57200 57400 57600 57800 1 10F [m Jy ] J-band (REM and SMARTS) 55800 56000 56200 56400 56600 56800 57000 57200 57400 57600 57800 0 1 2 3 4 F [J y] POLAMI 229 GHz 55800 56000 56200 56400 56600 56800 57000 57200 57400 57600 57800 0 2 4 6 8 F [J y] Metsahovi 37 GHz 55800 56000 56200 56400 56600 56800 57000 57200 57400 57600 57800 MJD 0 2 4 6 8 F [J y] OVRO 15 GHz Fig. 1. Multiwavelength light curve of PKS 1510–089 between 2012 and 2017. From top to bottom: Fermi-LAT flux above 1 GeV; Swift-XRT flux 2-10 keV; U band flux from UVOT; KVA, SMARTS and MAPCAT optical flux in R-band; IR flux from REM and SMARTS in J band; radio 229 GHz flux measured by POLAMI; radio 37 GHz flux measured by Metsähovi; 15 GHz flux measured by OVRO. The red points show the observations within 12 h (or 3 days for the radio measurements) when MAGIC data have been taken during the time that Fermi-LAT flux is above 3 × 10−8cm−2s−1, while the blue points are observations in time bins with Fermi-LAT flux below this flux value (i.e. the low-state sample). IR, optical and UV data have been dereddened using Schlafly & Finkbeiner (2011). is ∼ 1.5 × 10−8cm−2s−1. We also include in the low-state sam- ple nights for which the Fermi-LAT flux did not reach TS of 4. The average 95% C.L. flux upper limit on those nights is also ∼ 3 × 10−8cm−2s−1. Article number, page 4 of 13 MAGIC Collaboration: V. A. Acciari et al.: Detection of VHE gamma rays from PKS 1510–089 during low state. 8.5− 8− 7.5− 7− 6.5− 6− 5.5− )-1s-2log(Flux>1GeV / cm 1 10 210 # of d ay s or m on th s Fig. 2. Distribution of flux above 1 GeV measured with Fermi-LAT in 30-day (thick line) or 1-day (thin line) bins. The vertical dashed line shows the value of the cut separating the low state. For the low-state spectrum we combine individual one day integration windows selected with flux fulfilling Eq. 1. The spec- trum, calculated above 100 MeV, is well described as a power law with spectral index 2.56 ± 0.04, with a TS= 1656.0. The possible curvature in the spectrum is investigated by fitting the spectrum with a logparabola which yields a TS= 1655.42 and negligible curvature (β = 0.06 ± 0.04). Therefore, no hint of spectral curvature was found during the low-state periods con- sidered in this analysis. The selection of Fermi-LAT observation days according to the flux> 1 GeV can bias the reconstructed av- erage spectrum in this energy range. To investigate such possible bias for the selected low-emission time periods we also calculate the spectral index in the energy range 0.1–1 GeV (not affected by the data selection) and obtain 2.41 ± 0.06. Moreover, the Fermi- LAT spectral energy points above 1 GeV are ∼ 25% lower than the extrapolation of the spectrum below 1 GeV, suggesting that indeed there is an up to 25% underestimation bias in the obtained Fermi-LAT flux above 1 GeV. 2.2. MAGIC MAGIC is a system of two imaging atmospheric Cherenkov telescopes. The telescopes are located in the Canary Islands, on La Palma (28.7◦N, 17.9◦W), at a height of 2200 m above sea level (Aleksić et al. 2016a). The large mirror dish diame- ter of 17 m, resulting in low energy threshold, makes it a per- fect instrument for studies of soft-spectrum sources such as FSRQs. As PKS 1510–089 is a southern source, only observable at zenith angle > 38◦, the corresponding trigger threshold would be & 90 GeV for a Crab nebula-like spectrum (Aleksić et al. 2016b), about 1.7 times larger than for the low zenith observa- tions. About 70% of the data of PKS 1510–089 was taken at the culmination, with zenith angle < 40◦. Moreover, PKS 1510– 089 is intrinsingly soft; the analysis energy threshold is only ∼ 80 GeV for a source with a spectral index of ∼ −3.3. Note also that the energy threshold of Cherenkov telescopes is not a sharp one and the unfolding procedure allows us to recon- struct the source spectrum slightly below the nominal value of the threshold. Between 2012 and 2017 the MAGIC telescopes observed PKS 1510–089 during 151 nights, out of which 115 passed at least partially the data quality selection cuts. We then select ]2 [ deg2θ 0 0.1 0.2 ev en ts N 0 2000 4000 6000 8000 10000 12000 14000 Time = 75.39 h 87.5± = 23031.2 off = 24726; NonN 180.0± = 1694.8 exN σSignificance (Li&Ma) = 9.53 Fig. 3. Distribution of θ2, the squared angular distance between the re- constructed arrival direction of individual events and the nominal source position (points) or background estimation region (gray filled area) for MAGIC observations of PKS 1510–089. Dashed line shows the value of the θ2 up to which the significance of the detection (see the inset text) is calculated. the nights corresponding to the Fermi-LAT periods fulfilling the Eq. 1 condition. Such procedure results in low-state data stacked from 76 nights, amounting to a total observation time of 75 hrs. The cut on the flux > 1 GeV excludes the MAGIC data reporting the detections of the two flares observed in 2015 (Ahnen et al. 2017a) and 2016 (Zacharias et al. 2017), as well as most of the data used for the detection during the high state of 2012 (Aleksić et al. 2014). The data were analyzed using MARS, the standard analysis package of MAGIC (Zanin et al. 2013; Aleksić et al. 2016b). Due to evolving telescope perfor- mance the data have been divided into 6 analysis periods. Within each analysis period proper Monte Carlo simulations are used for the analysis. At the last stage the analysis results from all the periods are merged together. This low-state data set shows a gamma-ray excess with a significance of 9.5σ (see Fig. 3). 2.3. Swift-XRT Since 2006, the source is monitored in the X-ray band by the XRT instrument on board the Neil Gehrels Swift Observatory (XRT, Burrows et al. 2004). In total, 243 raw images are pub- licly available in SWIFTXRLOG (Swift-XRT Instrument, Log)5. From those we selected 17 images based on the simultaneity to the GeV low-flux state and contemporaneousness with MAGIC observations. The standard Swift-XRT analysis6 is described in detail in Evans et al. (2009). The data are processed following the procedure described by Fallah Ramazani et al. (2017), as- suming a fixed equivalent Galactic hydrogen column density nH = 6.89×1020 cm−2 reported by Kalberla et al. (2005). We de- fined the source region as a circle of 20 pixels (∼ 47”) centered on the source, and a background region as a ring also centered on the source with inner and outer radii of 40 (∼94”) and 80 pixels (∼188”), respectively. In order to calculate the average low-state X-ray spectrum of PKS 1510–089 we have combined all selected individual Swift- XRT pointings (see the blue points in Fig. 1), adding up to a total exposure of 30 ks. The 2–10 keV flux measured during 5 https://heasarc.gsfc.nasa.gov/W3Browse/swift/ swiftxrlog.html 6 http://www.swift.ac.uk/analysis/xrt/ Article number, page 5 of 13 https://heasarc.gsfc.nasa.gov/W3Browse/swift/swiftxrlog.html https://heasarc.gsfc.nasa.gov/W3Browse/swift/swiftxrlog.html http://www.swift.ac.uk/analysis/xrt/ A&A proofs: manuscript no. pks1510-low those 17 pointings shows clear variability. Fitting the flux with a constant value yields χ2/Ndof = 84/16; however, the ampli- tude of the variability is moderate (the RMS of the points is about ∼ 30% of the mean flux). The average spectrum is well fitted (χ2/Ndof = 187.7/214) with a power law with an index of 1.382 ± 0.020 and F2−10 keV = 8.14+0.25 −0.19 × 10−12 erg cm−2s−1. The spectral index does not show significant variability (fit to constant yields χ2/Ndof=31.19/16, which translates to a chance probability of ∼ 1.2%). A harder-when-brighter trend is only hinted, with a Pearson’s rank coefficient for a linear correlation between flux and spectral index of 0.81 (2-10 keV) and 0.74 (0.3- 10 keV). 2.4. Optical observations PKS 1510–089 is regularly monitored as part of the Tuorla blazar monitoring program7 in the R band using a 35 cm Celestron tele- scope attached to the KVA (Kunglinga Vetenskapsakademi) tele- scope located at La Palma. The monitoring covers the period of 2012-2017 and the observations were mostly contemporaneous with the MAGIC observations of the source. The data analysis was performed with the semi-automatic pipeline using the stan- dard analysis procedures (Nilsson et al. in prep). The differen- tial photometry was performed using the comparison star mag- nitudes from Villata et al. (1997). Calar Alto data were acquired as part of the MAPCAT (Mon- itoring AGN with Polarimetry at the Calar Alto 2.2m Tele- scope) project8, see Agudo et al. (2012). The MAPCAT data pre- sented here were reduced following the procedure explained in Jorstad et al. (2010). Additionally, we used the publicly available data in the R band from the Small and Moderate Aperture Research Tele- scope System (SMARTS) instrument located at Cerro Tololo In- teramerican Observatory (CTIO) in Chile (Bonning et al. 2012), processed as described in Ahnen et al. (2017a). The KVA, MAPCAT and SMARTS R-band data have been corrected for Galactic extinction using AR = 0.217 (Schlafly & Finkbeiner 2011). In the optical range, PKS 1510– 089 shows mostly low emission throughout 2012–2014 and dur- ing 2017. Strong flares are seen in 2015 and 2016 at the times of high GeV state. The individual measurements performed in the optical range have very small statistical uncertainties, well be- low the variability observed during the selected low-state nights. Therefore, for the modelling, we use the average optical flux from 53 nights of observations (47 nights of observations with KVA, 3 with MAPCAT and 13 with SMARTS). We take as the uncertainty the RMS spread of the measurements. By applying such a procedure, we obtain that the mean optical flux during the low state is 1.55 ± 0.57 mJy. In band B we combine Swift- UVOT data (see the next section) with the SMARTS data, bring- ing the total number of observing nights to 20 and average flux of 1.22 ± 0.46 mJy. 2.5. Swift-UVOT The Ultraviolet/Optical Telescope (UVOT, Poole et al. 2008) is an instrument on board the Swift satellite operating in the 180– 600 nm wavelength range. The source counts were extracted from a circular region centered on the source with 5” radius, the background counts from an annular region centered on the source with inner and outer radius of 15” and 25” respectively. 7 http://users.utu.fi/kani/1m 8 http://www.iaa.es/~iagudo/_iagudo/MAPCAT.html The data calibration was done following Raiteri et al. (2010), where the effective wavelength, counts-to-flux conversion fac- tor, and Galactic extinction for each filter were calculated in an iterative procedure by taking into account the filter’s effec- tive area and the source’s spectral shape. The Galactic extinc- tion values derived from the re-calibration procedure are Av = 0.28 mag, Ab = 0.37 mag, Au = 0.44 mag, Aw1 = 0.63 mag, Am2 = 0.78 mag, Aw2 = 0.74 mag. The variability of the UV flux during the low-state nights is rather minor. The average flux of the quiescent state was derived in the same way as for optical data. The number of quasi-simultaneous UVOT obser- vations, contemporaneous to MAGIC observations during the Fermi-LAT low gamma-ray state is 9–13, depending on the filter. 2.6. IR We use infrared observations of PKS 1510–089 performed with the REM (Rapid Eye Mount, Zerbi et al. 2001; Covino et al. 2004) 60 cm diameter telescope located at La Silla, Chile. The observations were performed with J, H and Ks filters, with indi- vidual exposures ranging from 12 s to 30 s. After calibration us- ing neighbouring stars, the magnitudes were converted to fluxes using the zero magnitude fluxes from Mead et al. (1990). Addi- tionally, we used the publicly available data in the J and K bands from SMARTS (Bonning et al. 2012), processed as described in Ahnen et al. (2017a). Since the data were taken independently from MAGIC, a limited number of nights of MAGIC observations have quasisi- multaneous REM or SMARTS data. The data taken during the times classified as low state consist of 5 nights of REM data for H filter and 13 nights of REM or SMARTS data from J and K fil- ters. Moreover, one of the nights observed by SMARTS on MJD 57181 had a major IR flare where the flux increased by a factor of ∼5–6 with respect to the average flux value of the rest of the selected data. We nevertheless apply the same procedure as for R-band KVA data, averaging the IR flux over those low-state ob- servations, neglecting, however, the night of the IR flare. We ob- tain FK = 7.3±2.7 mJy, FH = 4.2±2.4 mJy, FJ = 2.3±1.0 mJy. Including the night of the IR flare in the sample would change the FJ and FK fluxes relatively mildly (∼ 30% increase), however it would increase the RMS considerably to a value comparable to the flux. 2.7. Radio We use radio monitoring observations of PKS 1510–089 performed by OVRO (15 GHz, Richards et al. 2011), Met- sähovi (37 GHz, Teräesranta et al. 1998) and POLAMI (86 GHz, 229 GHz). We also use CARMA data taken at 95 GHz between August 2012 and November 2014 (Ramakrishnan et al. 2016). POLAMI (Polarimetric Monitoring of AGN at Millime- tre Wavelengths)9 (Agudo et al. 2018a,b; Thum et al. 2018) is a long-term program to monitor the polarimetric properties (Stokes I, Q, U, and V) of a sample of ∼40 bright AGN at 3.5 and 1.3 millimeter wavelengths with the IRAM 30m Telescope near Granada, Spain. The program has been running since Octo- ber 2006, and it currently has a time sampling of ∼2 weeks. The XPOL polarimetric observing setup has been routinely used as described in Thum et al. (2008) since the start of the program. The reduction and calibration of the POLAMI data presented here are described in detail in Agudo et al. (2010, 2014, 2018a). 9 http://polami.iaa.es Article number, page 6 of 13 http://users.utu.fi/kani/1m http://www.iaa.es/~iagudo/_iagudo/MAPCAT.html http://polami.iaa.es MAGIC Collaboration: V. A. Acciari et al.: Detection of VHE gamma rays from PKS 1510–089 during low state. The 37 GHz observations were made with the 13.7 m di- ameter Aalto University Metsähovi radio telescope10, which is a radome enclosed Cassegrain type antenna situated in Finland. The measurements were made with a 1 GHz-band dual beam re- ceiver centered at 36.8 GHz. The HEMPT (High Electron Mobil- ity Pseudomorphic Transistor) front end operates at room tem- perature. The observations are Dicke switched ON–ON obser- vations, alternating between the source and the sky in each feed horn. A typical integration time to obtain one flux density data point is between 1200 and 1800 s. The detection limit of the telescope at 37 GHz is of the order of 0.2 Jy under optimal con- ditions. Data points with a signal-to-noise ratio < 4 are consid- ered non-detections. The flux density calibration is set by ob- servations of the HII region DR 21. The sources NGC 7027, 3C 274 and 3C 84 are used as secondary calibrators. A de- tailed description of the data reduction and analysis is given in Teräesranta et al. (1998). The error estimate in the flux density includes the contribution from the measurement RMS and the uncertainty of the absolute calibration. The Owens Valley Radio Observatory (OVRO) 40-Meter Telescope uses off-axis dual-beam optics and a cryogenic re- ceiver with 3 GHz bandwidth centered at 15 GHz. Atmospheric and ground contributions as well as gain fluctuations are re- moved with the double switching technique (Readhead et al. 1989) where the observations are conducted in an ON-ON fash- ion so that one of the beams is always pointed on the source. Until May 2014 the two beams were rapidly alternated using a Dicke switch. Since then a new pseudo-correlation receiver re- placed the old one, and a 180◦ phase switch is used. Relative calibration is obtained with a temperature-stable noise diode to compensate for gain drifts. The primary flux density calibrator for those observations was 3C 286 with an assumed value of 3.44 Jy(Baars et al. 1977), while DR21 is used as a secondary calibrator source. Details of the observation and data reduction procedures are given in Richards et al. (2011). The radio flux at all frequencies shows slow variability, not simultaneous with the flares observed at higher energies. In or- der to obtain the average emission during the low gamma-ray state we apply the same procedure as for the R-band flux; how- ever we apply a larger margin in time, using the data within ±3 days from the MAGIC observations during low Fermi-LAT flux. We obtain F15 GHz = 4.4± 1.2 Jy (average over 22 observations), F37 GHz = 3.9±1.1 Jy (59 observations), F86 GHz = 3.14±0.86 Jy (6 observations), F95 GHz = 2.16 ± 0.13 Jy (9 observations), F229 GHz = 1.76 ± 0.42 Jy (4 observations). 3. Low gamma-ray state of PKS 1510–089 The low-state spectrum of PKS 1510–089 observed by the MAGIC telescopes was reconstructed between 63 and 430 GeV and is shown in Fig. 4. The observed spectrum can be de- scribed by a power law: dN/dE = (4.66 ± 0.59stat) × 10−11(E/175 GeV)−3.97±0.23statcm−2s−1TeV−1. Correcting for the absorption due to the interaction with the extragalactic back- ground light (EBL) according to Domínguez et al. (2011), we obtain the following intrinsic spectrum: dN/dE = (7.9±1.1stat)× 10−11(E/175 GeV)−3.26±0.30statcm−2s−1TeV−1. Since the excess to residual background ratio is of the order of 6%, the systematic uncertainty on the flux normalization (without the effect of the energy scale) is ±20%, larger than for typical MAGIC obser- vations, (Aleksić et al. 2016b). Also the systematic uncertainty on the spectral slope is increased by the low excess to residual 10 http://metsahovi.aalto.fi/en/ 70 80 90100 200 300 400 500 600 700 E [GeV] 14−10 13−10 12−10 11−10 10−10 ] -1 s -2 d N /d E [T eV c m 2 E MAGIC, 2012-2017 low-Fermi-LAT H.E.S.S., Mar 2009 MAGIC, Feb-Mar 2012 MAGIC, May 2015 flare MAGIC, May 2016 flare Fig. 4. Spectral Energy Distribution (SED) of PKS 1510–089 during the low state (red filled points and shaded magenta region) compared to historical measurements (open symbols): high state in March 2009 (gray stars, Abramowski et al. 2013), high state in February-March 2012 (black diamonds, Aleksić et al. 2014), flare in May 2015 (blue crosses, Ahnen et al. 2017a) and flare in May 2016 (magenta squares, (Zacharias et al. 2017)). The spectra are not deabsorbed from the EBL extinction. background ratio and following the prescription of Section 5.1 in Aleksić et al. (2016b) can be estimated as ±0.4. The uncer- tainty of the energy scale is ±15%. Comparing with previous measurements, the high-state detection in 2012 (Aleksić et al. 2014) gives ∼ 1.7 times larger flux than the low state stud- ied here. On the other hand, the most luminous flare observed from PKS 1510–089 in May 2016 gives flux ∼40–80 times higher than the low state (for the MAGIC and H.E.S.S. obser- vation window respectively, see Zacharias et al. 2017). Interest- ingly, the intrinsic spectral index of −3.26± 0.30stat is consistent within the uncertainties with the one obtained during the high state in the 2012 (−2.5 ± 0.6stat, Aleksić et al. 2014), the 2015 flare (−3.17 ± 0.80stat, Ahnen et al. 2017a) and the 2016 flare (−2.9 ± 0.2stat, −3.37 ± 0.09stat, Zacharias et al. 2017). As reported in Section 2, the IR to UV low-state data show variability at the level of ∼ 40%. We search for possible variabil- ity in the MAGIC data taken during the defined low state by com- puting light curves using different binnings (see Fig. 5). Both the daily (χ2/Ndof = 51.9/74) and yearly (χ2/Ndof = 3.08/5) light curves do not show any evidence of variability when fitted with a constant flux model. The gamma-ray flux of PKS 1510–089 dur- ing the low state is, however, too weak for probing variability with a similar relative amplitude at GeV energies with MAGIC or Fermi-LAT as observed in IR-UV. The average emission of the low state above 150 GeV is (4.27± 0.61stat)× 10−12 cm−2s−1, which is also below the all-time average of the H.E.S.S. obser- vations ((5.5 ± 0.4stat) × 10−12 cm−2s−1, Zacharias et al. 2016). 4. Modelling The multiwavelength SED constructed from the data selected ac- cording to the low flux above 1 GeV, taken between 2012 and 2017 is shown in Fig. 6. We model the broad-band emission us- ing an External Compton scenario (see, e.g., Sikora et al. 1994; Ghisellini et al. 2010) in which the gamma-ray emission is pro- duced due to inverse Compton scattering of a radiation field ex- ternal to the jet by electrons located in an emission region inside the jet. We use a particular scenario applied already to model a Article number, page 7 of 13 http://metsahovi.aalto.fi/en/ A&A proofs: manuscript no. pks1510-low MJD 56000 56500 57000 57500 58000 ) fo r E > 1 50 G eV -1 s -2 F lu x (c m 15− 10− 5− 0 5 10 15 20 25 12−10× Fig. 5. Light curve of PKS 1510–089 obtained with MAGIC observations during the low state in daily (red thin lines) and yearly (black thick lines) binning. For clarity three nights with short exposure (resulting in flux estimation uncertainty > 15×10−12cm−2s−1) are omitted from the plot. "close" "far" Fig. 6. Multiwavelength SED of PKS 1510–089 obtained from the data contemporaneous to MAGIC observations performed during Fermi-LAT low state (red points). The gray band shows the Swift-BAT 105 months average spectrum (Oh et al. 2018). Gray dot markers show the historical data from SSDC (www.asdc.asi.it). IR optical and UV data have been dereddened, MAGIC data have been corrected for the absorption by the EBL according to Domínguez et al. (2011) model. Observed MAGIC spectral points are shown in cyan. The green short-long-dashed curve shows the synchrotron component, and orange dot-dashed the EC component. SSC component is shown in cyan dotted line. The long dashed and short dashed black lines show the dust torus and accretion disk emission respectively. The solid blue line shows the total emission (including absorption in EBL). The left panel is for the “close” model, the right panel for the “far” model (see the text). high state and a flare from PKS 1510–089 (Aleksić et al. 2014; Ahnen et al. 2017a), with the external photon field being the ac- cretion disk radiation reflected by the broad line region (BLR) and dust torus (DT). We apply the same BLR and DT parameters as in Ahnen et al. (2017a), namely a radius of RBLR = 2.6 × 1017 cm and RIR = 6.5 × 1018 cm respectively. BLR and DT reflect fBLR = 0.1 and fDT = 0.6 respectively (the so-called covering factor) part of the accretion disk radiation, Ldisk = 6.7× 1045 erg s−1. The DT temperature is set to 1000 K. The emission region is located at the distance r from the base of the jet and has a radius of R. As in the model employed in Ahnen et al. (2017a), jet plasma flows through the emission region. The lack of strong variability of the low-state emission and the fact that it reaches sub-TeV energies suggests that the emission region should be beyond the BLR. At such distances the cross section of the jet is large, making difficult to explain any short-term variability, but the absorption on the BLR radiation is negligible. We con- sider two scenarios for the location of the emission region: the “close” scenario with r = 7 × 1017 cm and “far” scenario with r = 3 × 1018 cm. In both cases, the dominating radiation field comes from the DT. Such distances of the emission region have been applied for modelling of the 2015 flare (Ahnen et al. 2017a) and 2012 high state (Aleksić et al. 2014) respectively. The size of the emission region R = 2×1016 cm (for the “close” scenario) and R = 3 × 1017 cm (for the “far” scenario) is on the order of the cross section of the jet at the distance r. Although it is not Article number, page 8 of 13 www.asdc.asi.it MAGIC Collaboration: V. A. Acciari et al.: Detection of VHE gamma rays from PKS 1510–089 during low state. a dominant emission component, the model also calculates the synchrotron-self-Compton emission of the source. The model parameters for both scenarios are summarized and compared with earlier modelling in Table 1. However, the sets of parameters are not unique solutions for describing the low-state SED, as a certain degree of parameter degeneracy oc- curs in these kind of models (see e.g. synchrotron-self-Compton (SSC) model parameters degeneracy discussed in Ahnen et al. 2017b). The data are compared with the model in Fig. 6. Both sce- narios can reproduce relatively well the IC peak. The gamma- ray data of MAGIC and Fermi-LAT are well explained as the high-energy part of the EC component, with the exception of the two highest energy Fermi-LAT points which are > 1 GeV and hence are probably underestimated by the data selection proce- dure (see Section 2.1). Swift-XRT and historical Swift-BAT data follow the rising part of the EC component (with a small contri- bution of the SSC process in the soft X-ray range for the “close” scenario). The UV data form a bump that can be well explained by the direct accretion disk radiation included in the model. In the IR range, the model curve underestimates the data points, especially in the case of the “far” model. Among the quiescent data selected, the IR data show the highest variability. The higher IR variability might come from a separate region, not associated with the GeV gamma-ray emission region. In such a case, the IR emission associated with the low gamma-ray state would likely be at the level closer to the low edge of the observed spread in IR fluxes (reflected in the quoted uncertainty bar in Fig. 6). The “far” model can reasonably reproduce the radio obser- vations, while the “close” model underestimates the data due to strong synchrotron-self-absorption effects given by the compact- ness of the emission region. This is not surprising since the radio core observed at 15 GHz is estimated to be located at the dis- tance of 17.7 pc from the base of the jet (Pushkarev et al. 2012). Using the typical scaling of the core distance being inversely proportional to the frequency, we obtain that for the highest ra- dio point at 229 GHz its corresponding radio core should be lo- cated at ∼ 1 pc. Therefore, most of the radio emission should be produced at or beyond the region considered in the “far” sce- nario. However, the magnetic field considered in the “far” sce- nario, B = 0.05 G, is an order of magnitude smaller that the magnetic field estimated from the radio observations at r = 1 pc of 0.73 G (Pushkarev et al. 2012). Larger values would result in a much smaller Compton dominance than observed in the broad- band SED. It is curious that an optical/GeV high state, a days-long flare and the low state can all be roughly described (except of the IR data) in the framework of the same External Compton scenario without a major change of the model parameters. This suggests a common origin of the gamma-ray emission of PKS 1510–089 in different activity states, with the observed differences caused by changes in the content of the plasma flowing through the emis- sion region11. We note, however, that the model used here is rather simple. It is natural to assume that the low-state, broad- band emission is an integral of the emission in a range of dis- tances from the base, with the varying conditions (such as B field) along the jet, rather than originating in a single homoge- nous region (see e.g. Potter & Cotter 2013). 11 Note that a fast flare observed in 2016 from PKS 1510–089 (Zacharias et al. 2017) might have nevertheless a different origin. 1−10 1 10 210 E[GeV] 12−10 11−10 10−10 ] -1 s -2 S E D [T eV c m Fermi-LAT Out =1.00 RemR MAGIC Out =0.82 RemR MAGIC syst. Out =0.74 RemR Fig. 7. Gamma-ray spectrum of PKS 1510–089 during low state mea- sured by Fermi-LAT (squares) and MAGIC (filled circles). 68% con- fidence band of the extrapolation of the Fermi-LAT spectrum to sub- TeV energies is shown with gray shaded region. The extrapolation in the MAGIC energy range assuming absorption in BLR following (Böttcher & Els 2016) for the emission region located at the distance of 1, 0.82 and 0.74 of the outer radius of the broad line region is shown with black solid, blue dotted and green dashed lines respectively. Empty circles show the effect of the systematic uncertainties on the MAGIC spectrum. 4.1. Limits on the absorption of sub-TeV photons in BLR In the case the emission region is located inside, or close to, the BLR the gamma-ray spectrum should carry an imprint of the absorption feature on the BLR photons (Donea & Protheroe 2003). Presence or lack of such an absorption can be therefore used to constrain the location of the emission region. We use the emission-model-independent approach of Cerruti et al. (2017) to put such constraints of the location of the low-state emission region of PKS 1510–089. We first make a power-law fit to the Fermi-LAT spectrum of PKS 1510–089 in the energy range of 0.1–1 GeV, which is unbiased by the data selection. Next, we extrapolate the fit to the energy range observed by the MAGIC telescopes and apply an absorption by a factor of exp(−τ), where τ is the so-called optical depth (see Fig. 7). We compare those extrapolations with the reconstructed MAGIC spectrum taking into account the statistical uncertain- ties as well as the systematic uncertainty both in the energy scale and in the flux normalization. The systematic uncertain- ties of Fermi-LAT are negligible in those calculations. Due to source intrinsic effects (e.g. intrinsic break or cut-off in the ac- celerated electron spectrum, Klein-Nishina effect) the sub-TeV spectrum might be below the GeV extrapolation. Hence no lower limit on the absorption can be derived in a model-independent way. However it is natural to assume that there is no upturn- ing in the photon spectrum, therefore we can place an upper limit on the maximum absorption in the BLR. To estimate such a limit we perform a toy Monte Carlo study in which we vary 10000 times the extrapolated Fermi-LAT flux and the measured MAGIC flux (corrected for the EBL absorption) according to their uncertainties. Next, at each investigated energy we calcu- late a histogram of τ(E) = ln(Fext(E)/Fobs(E)), where Fext(E) and Fobs(E) are the randomized extrapolated and randomized measured flux respectively. The limit on τ is obtained as a 95% quantile of such a distribution. We include the systematic un- certainties of MAGIC by shifting its energy scale and normal- Article number, page 9 of 13 A&A proofs: manuscript no. pks1510-low γmin γb γmax n1 n2 B K δ Γ r R Low State (close) 2.5 130 3 × 105 1.9 3.5 0.35 3 × 104 25 20 7.0 × 1017 2.0 × 1016 Low State (far) 2 300 3 × 105 1.9 3.7 0.05 80 25 20 3.0 × 1018 3.0 × 1017 2012 3 900 6.5 × 104 1.9 3.85 0.12 20 20 20 3.1 × 1018 3.0 × 1017 2015, Period A 1 150 & 800 4 × 104 1 & 2 3.7 0.23 3.0 × 104 25 20 6.0 × 1017 2.8 × 1016 2015, Period B 1 150 & 500 3 × 104 1 & 2 3.7 0.34 2.6 × 104 25 20 6.0 × 1017 2.8 × 1016 Table 1. Input model parameters for the EC models of PKS 1510–089 emission for the low state in “close” and “far” scenario. For comparison, model parameters obtained from the 2012 high state (Aleksić et al. 2014) and 2015 flare Ahnen et al. (2017a) are also quoted. The individual columns are minimum, break and maximum electron Lorentz factor (γmin, γb, γmax respectively), slope of the electron energy distribution below and above γb (n1 and n2 respectively), magnetic field in G (B), normalization of the electron distribution in units of cm−3 (K), Doppler, Lorentz factor, distance and radius of the emission region (δ, Γ, r (in cm), R (in cm) respectively). ization according to their systematic limits (see empty circles in Fig. 7) and taking the least constraining one. We obtain that τ(110 GeV)<1.4, τ(180 GeV)<1.7, τ(290 GeV)<2.3. Applying a model of absorption by the BLR those limits on the optical depth can be converted to lower limits on the loca- tion of the emission region, r. We test the above procedure using two BLR models for PKS 1510–089. First, we use the optical depth calculations of Böttcher & Els (2016) assuming that 10% of the accretion disk radiation is reprocessed in the BLR. Note that Böttcher & Els (2016) assumes that a homogeneous BLR in PKS 1510–089 spans between 6.9 × 1017 cm and 8.4 × 1017 cm, and reflects 10% of the disk luminosity LD = 1046 erg s−1. We obtain that the above limits result in r > 6.3 × 1017 cm (i.e. above 0.74 of the outer radius of the BLR). As an additional check we calculate the optical depths using a code adapted from Sitarek & Bednarek (2008) with the line intensities and BLR ge- ometry of Liu & Bai (2006). We use the same radius and lumi- nosity of the BLR as in Section 4. We apply however the same ratio of the outer to inner radius of BLR as in Böttcher & Els (2016), resulting in the BLR spanning from 2.34 × 1017 cm to 2.86 × 1017 cm. For such a BLR model, the obtained above lim- its on the optical depth force us to place the emission region farther than 3.2 × 1017 cm (i.e. beyond 1.1 of the outer radius of the BLR. It should be noted that the method has a number of simpli- fications and underlying assumptions. The emission region is assumed to be relatively small compared to its distance from the black hole. This is not necessarily true, in particular for the low state emission that can be generated in a more ex- tended region (the broad band emission modelling presented in the previous section further supports such a hypothesis). Sec- ond, if the gamma-ray emission is not produced by a single pro- cess, the spectrum can have a complicated (including convex) shape. Note that for another FSRQ, B0218+357 the gamma-ray emission was explained as combination of SSC and EC process (Ahnen et al. 2016). Third, the optical depth is dependent on the assumed geometry of the BLR. For example, the size of the BLR derived in (Böttcher & Els 2016) is a factor of ∼ 3 larger than the one of Ghisellini & Tavecchio (2009). In addition the difference between the spherical and disk like geometry can easily change the optical depths by a factor of a few (Abolmasov & Poutanen 2017). Finally, the radial stratification of the BLR and the total fraction of the light reflected in the BLR introduce further uncer- tainties. 5. Conclusions We have performed the analysis of MAGIC data searching for a possible low state of VHE gamma-ray emission from PKS 1510– 089. Selecting the data taken during periods when the Fermi- LAT flux above 1 GeV was below 3× 10−8cm−2s−1 we have col- lected 75 hrs of MAGIC data on 76 individual nights, resulting in a significant detection of the low state of VHE gamma-ray emis- sion. The measured flux is ∼ 0.6 of the flux of the source mea- sured during high optical and GeV state (Aleksić et al. 2014) in the beginning of 2012 and∼ 0.75 of the lowest previously known flux from this source (average over all the H.E.S.S. observations, Zacharias et al. 2016). Nevertheless, the spectral shape is consis- tent with the previous measurements, despite a factor of 80 dif- ference to the flux during the strongest flare observed so far from this source. This makes PKS 1510–089 the first FSRQ to be de- tected in a persistent low state with no hints of yearly variations in the observed flux. Future observations with the Cherenkov Telescope Array should be able to probe if the low-state flux is also stable on shorter time scales (Acharya et al. 2013, 2017). Previous VHE gamma-ray observations of FSRQs in flar- ing states suggested that the emission region during such states should be located beyond the BLR and that the emission is mostly compatible with an EC scenario. The low-state broad- band emission of PKS 1510–089 from IR to VHE gamma-rays can be explained in the framework of an EC scenario, similarly to the previous VHE gamma-ray detections of the source. The presence of the sub-TeV gamma-ray emission also suggests that the emission region is located beyond the BLR, where the dom- inating seed radiation field for the EC process is the the dust torus. Comparison of the extrapolated Fermi-LAT spectrum and the MAGIC measured spectrum using two distinct BLR models allows us to put a limit on the location of the emission region be- yond the 0.74 of the outer radius of BLR. The emission scenario placing the dissipation region beyond the BLR is in line with the recent studies of Costamante et al. (2018) showing that most of the Fermi-LAT detected blazars (including PKS 1510–089) have GeV emission consistent with lack of BLR absorption. Acknowledgements. We would like to thank the Instituto de Astrofísica de Ca- narias for the excellent working conditions at the Observatorio del Roque de los Muchachos in La Palma. The financial support of the German BMBF and MPG, the Italian INFN and INAF, the Swiss National Fund SNF, the ERDF under the Spanish MINECO (FPA2015-69818-P, FPA2012-36668, FPA2015- 68378-P, FPA2015-69210-C6-2-R, FPA2015-69210-C6-4-R, FPA2015-69210- C6-6-R, AYA2015-71042-P, AYA2016-76012-C3-1-P, ESP2015-71662-C2-2- P, CSD2009-00064), and the Japanese JSPS and MEXT is gratefully ac- knowledged. This work was also supported by the Spanish Centro de Ex- celencia “Severo Ochoa” SEV-2012-0234 and SEV-2015-0548, and Unidad de Excelencia “María de Maeztu” MDM-2014-0369, by the Croatian Sci- ence Foundation (HrZZ) Project IP-2016-06-9782 and the University of Rijeka Project 13.12.1.3.02, by the DFG Collaborative Research Centers SFB823/C4 and SFB876/C3, the Polish National Research Centre grant UMO- 2016/22/M/ST9/00382 and by the Brazilian MCTIC, CNPq and FAPERJ. The Fermi LAT Collaboration acknowledges generous ongoing support from a num- ber of agencies and institutes that have supported both the development and the operation of the LAT as well as scientific data analysis. These include the Na- tional Aeronautics and Space Administration and the Department of Energy in the United States, the Commissariat à l’Energie Atomique and the Centre Na- Article number, page 10 of 13 MAGIC Collaboration: V. A. Acciari et al.: Detection of VHE gamma rays from PKS 1510–089 during low state. tional de la Recherche Scientifique / Institut National de Physique Nucléaire et de Physique des Particules in France, the Agenzia Spaziale Italiana and the Istituto Nazionale di Fisica Nucleare in Italy, the Ministry of Education, Culture, Sports, Science and Technology (MEXT), High Energy Accelerator Research Organiza- tion (KEK) and Japan Aerospace Exploration Agency (JAXA) in Japan, and the K. A. Wallenberg Foundation, the Swedish Research Council and the Swedish National Space Board in Sweden. Additional support for science analysis during the operations phase is gratefully acknowledged from the Istituto Nazionale di Astrofisica in Italy and the Centre National d’Études Spatiales in France. This work performed in part under DOE Contract DE-AC02-76SF00515. This paper has made use of up-to-date SMARTS optical/near-infrared light curves that are available at www.astro.yale.edu/smarts/glast/home.php. IA acknowl- edges support by a Ramón y Cajal grant of the Ministerio de Economía, In- dustria, y Competitividad (MINECO) of Spain. Acquisition and reduction of the POLAMI and MAPCAT data was supported in part by MINECO through grants AYA2010-14844, AYA2013-40825-P, and AYA2016-80889-P, and by the Re- gional Government of Andalucía through grant P09-FQM-4784. The POLAMI observations were carried out at the IRAM 30m Telescope. IRAM is supported by INSU/CNRS (France), MPG (Germany) and IGN (Spain). The MAPCAT observations were carried out at the German-Spanish Calar Alto Observatory, which is jointly operated by the Max-Plank-Institut für Astronomie and the Insti- tuto de Astrofísica de Andalucía-CSIC This research has made use of data from the OVRO 40-m monitoring program (Richards et al. 2011) which is supported in part by NASA grants NNX08AW31G, NNX11A043G, and NNX14AQ89G and NSF grants AST-0808050 and AST-1109911. This publication makes use of data obtained at the Metsähovi Radio Observatory, operated by Aalto Univer- sity, Finland. The authors would like to thank a number of people who provided comments to the manuscript: R. Angioni, N. MacDonald, M. Giroletti, R. Ca- puto, D. Thompson and the anonymous journal reviewer. We would also like to thank M. Böttcher for providing us the numerical values of the optical depths calculated in Böttcher & Els (2016). References Abdo, A. A., Ackermann, M., Agudo, I., et al. 2010, ApJ, 721, 1425 Abolmasov, P., & Poutanen, J. 2017, MNRAS, 464, 152 H.E.S.S. Collaboration, Abramowski, A., Acero, F., et al. 2013, A&A, 554, A107 Acero, F., Ackermann, M., Ajello, M. et al. 2015, ApJS, 218, 23 Acero, F., Ackermann, M., Ajello, M., et al. 2016, ApJS, 223, 26 Acharya, B. S., Actis, M., Aghajani, T., et al. 2013, Astroparticle Physics, 43, 3 Acharya, B. S., Agudo, I., Al Samarai, I., et al. 2017, arXiv:1709.07997 Agudo, I., Thum, C., Wiesemeyer, H., & Krichbaum, T. P. 2010, ApJS, 189, 1 Agudo, I., Molina, S. N., Gómez, J. L., et al. 2012, International Journal of Mod- ern Physics Conference Series, 8, 299 Agudo, I., Thum, C., Gómez, J. L., & Wiesemeyer, H. 2014, A&A, 566, A59 Agudo, I., Thum, C., Molina, S. N., et al. 2018, MNRAS, 474, 1427 Agudo, I., Thum, C., Ramakrishnan, V., et al. 2018, MNRAS, 473, 1850 Ahnen, M. L., Ansoldi, S., Antonelli, L. A., et al. 2016, A&A, 595, A98 Ahnen, M. L., Ansoldi, S., Antonelli, L. A., et al. 2017a, A&A, 603, A29 Ahnen, M. L., Ansoldi, S., Antonelli, L. A., et al. 2017b, A&A, 603, A31 Aleksić, J., Ansoldi, S., Antonelli, L. A., et al. 2014, A&A, 569, A46 Aleksić, J., Ansoldi, S., Antonelli, L. A., et al. 2016, Astroparticle Physics, 72, 61 Aleksić, J., Ansoldi, S., Antonelli, L. A., et al. 2016, Astroparticle Physics, 72, 76 Atwood, W. B., Abdo, A. A., Ackermann, M., et al. 2009, ApJ, 697, 1071 Baars, J. W. M., Genzel, R., Pauliny-Toth, I. I. K., & Witzel, A. 1977, A&A, 61, 99 Bonning, E., Urry, C. M., Bailyn, C., et al. 2012, ApJ, 756, 13 Böttcher, M., & Els, P. 2016, ApJ, 821, 102 Brown, A. M. 2013, MNRAS, 431, 824 Burrows, D. N., Hill, J. E., Nousek, J. A., et al. 2004, Proc. SPIE, 5165, 201 Cerruti, M., Lenain, J.-P., Prokoph, H., & for the H. E. S. S. Collaboration 2017, arXiv:1708.00658 Proc of 35th ICRC, Busan, Korea Costamante, L., Cutini, S., Tosti, G., Antolini, E., & Tramacere, A. 2018, MN- RAS, arXiv:1804.02408 Covino, S., Stefanon, M., Sciuto, G., et al. 2004, Proc. SPIE, 5492, 1613 Domínguez, A., Primack, J. R., Rosario, D. J., et al. 2011, MNRAS, 410, 2556 Donea, A.-C., & Protheroe, R. J. 2003, Astroparticle Physics, 18, 377 Evans, P. A., Beardmore, A. P., Page, K. L., et al. 2009, MNRAS, 397, 1177 Fallah Ramazani, V., Lindfors, E., & Nilsson, K. 2017, A&A, 608, A68 Ghisellini, G., & Tavecchio, F. 2009, MNRAS, 397, 985 Ghisellini, G., Tavecchio, F., Foschini, L., et al. 2010, MNRAS, 402, 497 Jorstad, S. G., Marscher, A. P., Larionov, V. M., et al. 2010, ApJ, 715, 362 Kalberla, P. M. W., Burton, W. B., Hartmann, D., et al. 2005, A&A, 440, 775 Liller, M. H., & Liller, W. 1975, ApJ, 199, L133 Liu, H. T., & Bai, J. M. 2006, ApJ, 653, 1089 Marscher, A. P., Jorstad, S. G., Larionov, V. M., et al. 2010, ApJ, 710, L126 Mattox J. R., Bertsch, D. L., Chiang, J., et al., 1996, ApJ, 461, 396 Mead, A. R. G., Ballard, K. R., Brand, P. W. J. L., et al. 1990, A&AS, 83, 183 Oh, K., Koss, M., Markwardt, C. B., et al. 2018, arXiv:1801.01882 Poole, T. S., Breeveld, A. A., Page, M. J., et al. 2008, MNRAS, 383, 627 Potter, W. J., & Cotter, G. 2013, MNRAS, 429, 1189 Prince, R., Majumdar, P., & Gupta, N. 2017, ApJ, 844, 62 Raiteri, C. M., Villata, M., Bruschini, L., et al. 2010, A&A, 524, A43 Pushkarev, A. B., Hovatta, T., Kovalev, Y. Y., et al. 2012, A&A, 545, A113 Ramakrishnan, V., Hovatta, T., Tornikoski, M., et al. 2016, MNRAS, 456, 171 Readhead, A. C. S., Lawrence, C. R., Myers, S. T., et al. 1989, ApJ, 346, 566 Richards, J. L., Max-Moerbeck, W., Pavlidou, V., et al. 2011, ApJS, 194, 29 Saito, S., Stawarz, Ł., Tanaka, Y. T., et al. 2013, ApJ, 766, L11 Schlafly, E. F. and Finkbeiner, D. P., 2011, ApJ, 737, 103S Sikora, M., Begelman, M. C., & Rees, M. J. 1994, ApJ, 421, 153 Sitarek, J., & Bednarek, W. 2008, MNRAS, 391, 624 Tanner, A. M., Bechtold, J., Walker, C. E., Black, J. H., & Cutri, R. M. 1996, AJ, 112, 62 Teraësranta, H., Tornikoski, M., Mujunen, A., et al. 1998, A&AS, 132, 305 Thum, C., Wiesemeyer, H., Paubert, G., Navarro, S., & Morris, D. 2008, PASP, 120, 777 Thum, C., Agudo, I., Molina, S. N., et al. 2018, MNRAS, 473, 2506 Villata M., Raiteri, C. M., Ghisellini, G. et al. 1997 A&AS 121, 119 Zacharias, M., Böttcher, M., Chakraborty, N., et al. 2016, arXiv:1611.02098 Zacharias, M., Sitarek, J., Dominis Prester, D., et al. 2017, arXiv:1708.00653 Zanin, R., Carmona, E., Sitarek, J., et al., 2013, Proc of 33rd ICRC, Rio de Janeiro, Brazil, Id. 773 Zerbi, R. M., Chincarini, G., Ghisellini, G., et al. 2001, Astronomische Nachrichten, 322, 275 1 Inst. de Astrofísica de Canarias, E-38200 La Laguna, and Universi- dad de La Laguna, Dpto. Astrofísica, E-38206 La Laguna, Tenerife, Spain 2 Università di Udine, and INFN Trieste, I-33100 Udine, Italy 3 National Institute for Astrophysics (INAF), I-00136 Rome, Italy 4 ETH Zurich, CH-8093 Zurich, Switzerland 5 Università di Padova and INFN, I-35131 Padova, Italy 6 Technische Universität Dortmund, D-44221 Dortmund, Germany 7 Croatian MAGIC Consortium: University of Rijeka, 51000 Rijeka, University of Split - FESB, 21000 Split, University of Zagreb - FER, 10000 Zagreb, University of Osijek, 31000 Osijek and Rud- jer Boskovic Institute, 10000 Zagreb, Croatia. 8 Saha Institute of Nuclear Physics, HBNI, 1/AF Bidhannagar, Salt Lake, Sector-1, Kolkata 700064, India 9 Max-Planck-Institut für Physik, D-80805 München, Germany 10 now at Centro Brasileiro de Pesquisas Físicas (CBPF), 22290-180 URCA, Rio de Janeiro (RJ), Brasil 11 Unidad de Partículas y Cosmología (UPARCOS), Universidad Com- plutense, E-28040 Madrid, Spain 12 University of Łódź, Department of Astrophysics, PL-90236 Łódź, Poland 13 Deutsches Elektronen-Synchrotron (DESY), D-15738 Zeuthen, Germany 14 Institut de Física d’Altes Energies (IFAE), The Barcelona Institute of Science and Technology (BIST), E-08193 Bellaterra (Barcelona), Spain 15 Università di Siena and INFN Pisa, I-53100 Siena, Italy 16 Università di Pisa, and INFN Pisa, I-56126 Pisa, Italy 17 Universität Würzburg, D-97074 Würzburg, Germany 18 Finnish MAGIC Consortium: Tuorla Observatory (Department of Physics and Astronomy) and Finnish Centre of Astronomy with ESO (FINCA), University of Turku, FIN-20014 Turku, Finland; As- tronomy Division, University of Oulu, FIN-90014 Oulu, Finland 19 Departament de Física, and CERES-IEEC, Universitat Autónoma de Barcelona, E-08193 Bellaterra, Spain 20 Universitat de Barcelona, ICCUB, IEEC-UB, E-08028 Barcelona, Spain 21 Japanese MAGIC Consortium: ICRR, The University of Tokyo, 277-8582 Chiba, Japan; Department of Physics, Kyoto University, Article number, page 11 of 13 www.astro.yale.edu/smarts/glast/home.php http://arxiv.org/abs/1709.07997 http://arxiv.org/abs/1708.00658 http://arxiv.org/abs/1804.02408 http://arxiv.org/abs/1801.01882 http://arxiv.org/abs/1611.02098 http://arxiv.org/abs/1708.00653 A&A proofs: manuscript no. pks1510-low 606-8502 Kyoto, Japan; Tokai University, 259-1292 Kanagawa, Japan; RIKEN, 351-0198 Saitama, Japan 22 Inst. for Nucl. Research and Nucl. Energy, Bulgarian Academy of Sciences, BG-1784 Sofia, Bulgaria 23 Humboldt University of Berlin, Institut für Physik D-12489 Berlin Germany 24 also at Dipartimento di Fisica, Università di Trieste, I-34127 Trieste, Italy 25 also at Port d’Informació Científica (PIC) E-08193 Bellaterra (Barcelona) Spain 26 also at INAF-Trieste and Dept. of Physics & Astronomy, University of Bologna 27 INAF, Osservatorio Astrofisico di Torino, via Osservatorio 20, I- 10025 Pino Torinese, Italy 28 Università dell’Insubria, Dipartimento di Scienza ed Alta Tecnolo- gia, Via Valleggio 11, I-22100, Como, Italy 29 INAF-Istituto Nazionale di Astrofisica, Osservatorio Astronomico di Brera, Via Bianchi 46, I-23807 Merate (LC), Italy 30 Owens Valley Radio Observatory, California Institute of Technol- ogy, Pasadena, CA 91125, USA 31 Tuorla Observatory, University of Turku, Väisäläntie 20, FI-21500 Piikkiö, Finland 32 Departamento de Astronomía, Universidad de Chile, Camino El Ob- servatorio 1515, Las Condes, Santiago, Chile 33 Aalto University Metsähovi Radio Observatory, Metsähovintie 114, 02540 Kylmälä, Finland 34 Aalto University Department of Electronics and Nanoengineering, P.O. BOX 15500, FI-00076 AALTO, Finland. 35 Tartu Observatory, Observatooriumi 1, 61602 Tõravere, Estonia 36 Instituto de Radio Astronomía Millimétrica, Avenida Divina Pas- tora, 7, Local 20, E–18012 Granada, Spain 37 Instituto de Astrofísica de Andalucía (CSIC), Apartado 3004, E– 18080 Granada, Spain 38 Max–Planck–Institut für Radioastronomie, Auf dem Hügel, 69, D– 53121, Bonn, Germany Article number, page 12 of 13 MAGIC Collaboration: V. A. Acciari et al.: Detection of VHE gamma rays from PKS 1510–089 during low state. Appendix A: Data selection tests As discussed in Section 2.1 the energy threshold of 1 GeV used for the selection was motivated by its proximity to the VHE energy range. However, the daily estimation of the flux above 1 GeV have large uncertainty which can in principle affect the data selection. To test this effect we have applied instead a cut to select nights with flux above 100 MeV measured by Fermi-LAT to be below 10−6cm−2s−1. Such a cut results in the same number of the MAGIC observations nights (76) selected for the low state analysis. The number of individual nights that would change the classification to the low state or to the high state with such a cut is 9 (corresponding to 12% of the low-state sample) each. We have tested the validity of the used here low-state data selection procedure by applying a cut at the flux above 100 MeV instead of 1 GeV. We therefore conclude, that the value of the Fermi-LAT analysis threshold do not have a large impact on the selection of nights used for this analysis. We have tested the effect of leav- ing the spectral index free in the light curve analysis, and found that this does not strongly affect the fraction of the data which results to be classified as low GeV state, following the definition described above. Fixing the spectral index to the average value of 2.36 (see the 3FGL catalog, Acero et al. 2015) would change the number of nights assigned to the low state and high state by . 1% each. Article number, page 13 of 13 1 Introduction 2 Data 2.1 Fermi-LAT 2.2 MAGIC 2.3 Swift-XRT 2.4 Optical observations 2.5 Swift-UVOT 2.6 IR 2.7 Radio 3 Low gamma-ray state of PKS1510–089 4 Modelling 4.1 Limits on the absorption of sub-TeV photons in BLR 5 Conclusions A Data selection tests