Draftversion January10,2012 PreprinttypesetusingLATEXstyleemulateapjv.5/2/11 PG 1553+113: FIVE YEARS OF OBSERVATIONS WITH MAGIC J. Aleksic´1, E. A. Alvarez2, L. A. Antonelli3, P. Antoranz4, M. Asensio2, M. Backes5, J. A. Barrio2, D. Bastieri6, J. Becerra Gonza´lez7,8, W. Bednarek9, A. Berdyugin10, K. Berger7,8, E. Bernardini11, A. Biland12, O. Blanch1, R. K. Bock13, A. Boller12, G. Bonnoli3, D. Borla Tridon13, I. Braun12, T. Bretz14,26, A. Can˜ellas15, E. Carmona13, A. Carosi3, P. Colin13, E. Colombo7, J. L. Contreras2, J. Cortina1, L. Cossio16, S. Covino3, F. Dazzi16,27, A. De Angelis16, G. De Caneva11, E. De Cea del Pozo17, B. De Lotto16, C. Delgado Mendez7,28, A. Diago Ortega7,8, M. Doert5, A. Dom´ınguez18, D. Dominis Prester19, D. Dorner12, M. Doro20, D. Elsaesser14, D. Ferenc19, M. V. Fonseca2, L. Font20, C. Fruck13, R. J. Garc´ıa Lo´pez7,8, M. Garczarczyk7, 2 D. Garrido20, G. Giavitto1, N. Godinovic´19, D. Hadasch17, D. Ha¨fner13, A. Herrero7,8, D. Hildebrand12, 1 D. Ho¨hne-Mo¨nch14, J. Hose13, D. Hrupec19, B. Huber12, T. Jogler13, H. Kellermann13, S. Klepser1, 0 T. Kra¨henbu¨hl12, J. Krause13, A. La Barbera3, D. Lelas19, E. Leonardo4, E. Lindfors10, S. Lombardi6, 2 M. Lo´pez2, A. Lo´pez-Oramas1, E. Lorenz12,13, M. Makariev21, G. Maneva21, N. Mankuzhiyil16, K. Mannheim14, L. Maraschi3, M. Mariotti6, M. Mart´ınez1, D. Mazin1,13, M. Meucci4, J. M. Miranda4, R. Mirzoyan13, n H. Miyamoto13, J. Moldo´n15, A. Moralejo1, P. Munar-Adrover15, D. Nieto2, K. Nilsson10,29, R. Orito13, a I. Oya2, D. Paneque13, R. Paoletti4, S. Pardo2, J. M. Paredes15, S. Partini4, M. Pasanen10, F. Pauss12, J M. A. Perez-Torres1, M. Persic16,22, L. Peruzzo6, M. Pilia23, J. Pochon7, F. Prada18, P. G. Prada Moroni24, 8 E. Prandini6,∗, I. Puljak19, I. Reichardt1, R. Reinthal10, W. Rhode5, M. Ribo´15, J. Rico25,1, S. Ru¨gamer14, A. Saggion6, K. Saito13, T. Y. Saito13, M. Salvati3, K. Satalecka2, V. Scalzotto6, V. Scapin2, C. Schultz6, ] T. Schweizer13, M. Shayduk13, S. N. Shore24, A. Sillanpa¨a¨10, J. Sitarek9, D. Sobczynska9, F. Spanier14, O S. Spiro3, V. Stamatescu1, A. Stamerra4, B. Steinke13, J. Storz14, N. Strah5, T. Suric´19, L. Takalo10, C H. Takami13, F. Tavecchio3,∗, P. Temnikov21, T. Terzic´19, D. Tescaro24, M. Teshima13, O. Tibolla14, D. F. Torres25,17, A. Treves23, M. Uellenbeck5, H. Vankov21, P. Vogler12, R. M. Wagner13, Q. Weitzel12, . h V. Zabalza15, F. Zandanel18, R. Zanin1, p and, - S. Buson6, D. Horan30, S. Larsson31,32,33, F. D’Ammando34 o 1 IFAE,EdificiCn.,CampusUAB,E-08193Bellaterra,Spain tr 2 UniversidadComplutense,E-28040Madrid,Spain s 3 INAFNationalInstitute forAstrophysics,I-00136Rome,Italy a 4 Universit`adiSiena,andINFNPisa,I-53100Siena,Italy [ 5 TechnischeUniversit¨atDortmund,D-44221Dortmund,Germany 6 Universit`adiPadovaandINFN,I-35131Padova, Italy 4 7 Inst. deAstrof´ısicadeCanarias,E-38200LaLaguna, Tenerife,Spain v 8 Depto. deAstrof´ısica,UniversidaddeLaLaguna, E-38206LaLaguna,Spain 9 UniversityofL ´od´z,PL-90236Lodz,Poland 4 10 TuorlaObservatory, UniversityofTurku,FI-21500Piikkio¨,Finland 6 11 Deutsches Elektronen-Synchrotron (DESY),D-15738Zeuthen, Germany 7 12 ETHZurich,CH-8093Switzerland 2 13 Max-Planck-Institutfu¨rPhysik,D-80805Mu¨nchen, Germany . 14 Universit¨atWu¨rzburg,D-97074Wu¨rzburg,Germany 1 15 UniversitatdeBarcelona(ICC/IEEC), E-08028Barcelona,Spain 0 16 Universit`adiUdine,andINFNTrieste,I-33100Udine,Italy 1 17 InstitutdeCi`enciesdel’Espai(IEEC-CSIC),E-08193Bellaterra,Spain 1 18 Inst. deAstrof´ısicadeAndaluc´ıa(CSIC), E-18080Granada, Spain : 19 CroatianMAGICConsortium,RudjerBoskovicInstitute, UniversityofRijekaandUniversityofSplit,HR-10000Zagreb,Croatia v 20 UniversitatAut`onomadeBarcelona,E-08193Bellaterra,Spain i 21 Inst. forNucl. ResearchandNucl. Energy,BG-1784Sofia,Bulgaria X 22 INAF/OsservatorioAstronomicoandINFN,I-34143Trieste,Italy r 23 Universit`adell’Insubria,Como,I-22100Como,Italy a 24 Universit`adiPisa,andINFNPisa,I-56126Pisa,Italy 25 ICREA,E-08010Barcelona,Spain 26 nowat: Ecolepolytechniquef´ed´eraledeLausanne(EPFL),Lausanne, Switzerland 27 supportedbyINFNPadova 28 nowat: CentrodeInvestigaciones Energ´eticas,Medioambientales yTecnolo´gicas (CIEMAT),Madrid,Spain 29 nowat: FinnishCentreforAstronomywithESO(FINCA),UniversityofTurku,Finland 30 LaboratoireLeprince-Ringuet,Ecolepolytechnique, CNRS/IN2P3,Palaiseau,France 31 DepartmentofPhysics,Stockholm University,AlbaNova,SE-10691Stockholm, Sweden 32 TheOskarKleinCenterforCosmoparticlePhysics,AlbaNova,SE-10691Stockholm, Sweden 33 Department ofAstronomy,Stockholm University,SE-10691Stockholm, Sweden 34 ConsorzioInteruniversitarioFisicaSpaziale(CIFS), v.leSettimioSevero63,I-10133Torino,Italy and ∗ correspondingauthors. E-mail:[email protected],[email protected]. Draft version January 10, 2012 ABSTRACT We present the results of five years (2005–2009) of MAGIC observations of the BL Lac object PG 1553+113at very high energies (VHEs, E > 100GeV). Power law fits of the individual years are compatible with a steady mean photon index Γ = 4.27 ± 0.14. In the last three years of data, the 2 Aleksi´c et al. flux level above 150GeV shows a clear variability (probability of constant flux < 0.001%). The flux variations are modest, lying in the range from 4% to 11% of the Crab Nebula flux. Simultaneous optical data also show only modest variability that seems to be correlated with VHE gamma ray variability. We also performed a temporal analysis of (all available) simultaneous Fermi/LAT data of PG 1553+113 above 1GeV, which reveals hints of variability in the 2008–2009 sample. Finally, we present a combination of the mean spectrum measured at very high energies with archival data available for other wavelengths. The mean spectral energy distribution can be modeled with a one– zone Synchrotron Self Compton (SSC) model, which gives the main physical parameters governing the VHE emission in the blazar jet. Subject headings: radiation mechanisms: non-thermal, gamma-rays: observations, BL Lacertae ob- jects: individual (PG 1553+113) 1. INTRODUCTION correcting the observed spectrum for the extragalactic background light absorption, cannot be harder than a The majority of extragalactic γ ray sources, both at fixed value given by theory. The technique, applied to GeV energies and above 100GeV, are blazars, radio- PG 1553+113, lead to an upper limit of z < 0.74. Re- loud active galactic nuclei with a relativistic jet point- cently, Prandini et al. (2010) extended this method us- ing towards the Earth. Their emission is dominated by ing the spectrum measured at lower energies as limiting the non-thermalcontinuum produced within the jet and slope for the original spectrum, obtaining z < 0.66 for boosted by relativistic effects (Urry & Padovani 1995). PG1553+113,at2σ level. Otherapproachesrequirethe The spectral energy distribution (SED) displays two absence of a pile up athigh energies. With this method, broad peaks, widely interpreted as due to synchrotron, Mazin & Goebel (2007) get z < 0.42. Hence, in the low frequency peak, and inverse Compton, high fre- case of PG 1553+113, the upper limits obtained with quencypeak,mechanism(althoughthehighenergypeak these methods are in the range of the limits set by op- could also be the result of hadronic processes, as pro- tical measurements. In this work we adopt the redshift posed in Mannheim 1993). Among blazars, BL Lac ob- z = 0.40. Such a large redshift is also supported by the jects are characterized by extremely weak emission lines absence of significant points at energies above 700GeV in their optical spectra, which often makes a measure- in the spectrum of the source. mentoftheir redshiftdifficult. The largemajorityofex- AtMeV-GeVenergies,PG1553+113wasnotdetected tragalactic sources detected above 100GeV are BL Lac byEGRET,butitiswellvisiblebytheLargeAreaTele- objects,inwhichthepeakofthesynchrotronbumpislo- scope (LAT) on-boardFermi, being detected with a sig- cated in the UV-X-ray bands and the high energy peak nificance above 10σ already in the first three months around 100GeV (these sources are often called high fre- ofobservations(Abdo et al.2009a). Abdo et al.(2010b) quency peaked BL Lacs, HBLs). PG 1553+113 is a show that the Fermi/LAT spectrum is surprisingly con- BL Lac discovered by Green et al. (1986). The large stant both in normalization and slope over ∼200 days. X-ray to radio flux ratio makes this source a typical Interestingly, a stability of the spectrum was also sug- HBL. Indeed, its synchrotron peak is located between gested by the H.E.S.S. and MAGIC observations, show- the UV and X-ray bands. Its optical spectrum is fea- ingrather marginalvariabilityduring 2005and2006ob- tureless, preventing the direct determination of the red- servations. The stability of the VHE γ ray emission is shift. Indirect methods based on the non detection of in contrast to the behavior commonly observed in other the characteristic lines and of the host galaxy provide TeV emitting BL Lacs, showing rather pronounced vari- lower limits, ranging from 0.09 to 0.78 (Sbarufatti et al. ations at all timescales. 2006; Sbarufatti, Treves, & Falomo 2005). The most re- After its discovery, PG 1553+113 was regularly ob- centestimate,basedontheLyalphaforestmethod,gives a z ∼ 0.40–0.45 (Danforth et al. 2010). served by MAGIC. In this paper, we present the analy- sis of the new data taken from 2007 to 2009, combined PG 1553+133 has been discovered as a VHE γ with previousobservations. The differentialandintegral ray emitter by H.E.S.S. (Aharonian et al. 2006a) and fluxesareanalyzed,incomparisonwithpartiallysimulta- MAGIC (Albert et al. 2007a), with a flux of approxi- neous measurements at other wavelengths, and the sta- mately 2% of that of the Crab Nebula above 200GeV. bility of the spectrum over this long period is studied. The spectrum appears extremely soft (photon index Γ ∼4), as expected by the absorption of VHE photons Finally,we combine allthe data availableandmodel the SED with a one–zone Synchrotron Self Compton (SSC) through interaction with the Extragalactic Background model, constraining the main physical parameters that Light (EBL) if the source is located at relatively large governthe VHE emission in the blazar jet. redshift (Stecker, de Jager, & Salamon 1992). The ab- sorption process, in fact, is a function of the energy of 2. MAGICOBSERVATIONSANDDATAANALYSIS the photon and of the distance it has traveled. Spectra with indices Γ∼4 havebeen observedin blazarslocated Since autumn 2009, MAGIC (Cortina et al. 2009) is a stereo system composed by two Imaging Atmospheric at redshifts above 0.2. Cherenkov Telescopes (IACTs) located on La Palma, VHE γ ray observations have been used as alter- ◦ ◦ Canary Islands, Spain (28.75 N, 17.89 W, 2240m asl). native method to constrain the distance of blazars. In this paper, we present only data collected before Aharonian et al. (2006b) proposed a way to set an up- the stereo upgrade, with a single telescope, MAGIC I per limit on the distance of blazars based on the as- (Baixeras et al. 2004), hereafter called MAGIC. The sumption that the VHE intrinsic spectrum, obtained by parabolic-shaped reflector, with a total mirror area of PG 1553+113: five years of observations with MAGIC 3 236m2, allows MAGIC to collect the Cherenkov light TABLE 2 and focus it onto a multi-pixel camera, composed of PG 1553+113SIGNAL 577 photo-multipliers. MAGIC camera and trigger are designed to record data not only during dark nights, Year Time Opt. PSF EnergyTh. Excesses Signif. but also under moderate light conditions (i.e. moder- [h] [mm] [GeV] [σ] ate moon, twilight). Due to its comparatively low trig- 2007 11.5 13 80 1400±242 5.8 ger energy threshold of ∼50GeV, MAGIC is well suited 2008 8.7 13 150 542±69 8.1 to perform multiwavelength observations together with 2009 8.5 14.8 160 212±52 4.2 instruments operating in the GeV range. aPG1553+113 signal study. Fromleftto right: year of observa- tion,effectivetimeofgoodqualitydatausedforthesignalanalysis, TABLE 1 opticalpointspreadfunction(PSF),energythresholdoftheanaly- PG 1553+113FINAL DATA SET sis,numberofexcesseventsobservedandsignificanceofthesignal. by bad weather (including calima, i.e. Saharan sand- Cycle Date Eff. Time Zd Rate DC [min] [◦] [Hz] [µA] dust in the atmosphere) that limited the final data set and resulted in an increased energy threshold. 23/03/2007 58 19-29 164 darknight All data analyzed here were taken in the false-source 19/04/2007 32 22-28 155 darknight tracking(wobble)mode(Fomin et al.1994),inwhichthe 20/04/2007 150 17-29 163 darknight III 21/04/2007 115 17-24 155 darknight telescope pointing was alternated every 20 minutes be- ◦ 22/04/2007 101 17-23 162 darknight tween two sky positions at 0.4 offset from the source. ◦ 23/04/2007 143 17-27 161 darknight The zenith angle of 2007 observations varied from 17 ◦ ◦ 24/04/2007 92 17-23 160 darknight to 30 , in 2008 it extended up to 36 , while in 2009 it 17/03/2008 58 17-19 150 darknight coveredthe range from 17◦ to 35◦. 18/03/2008 26 18-19 150 darknight The data were analyzed using the standard MAGIC 01/04/2008 43 20-26 167 darknight analysischain(Albert et al.2008a;Aliu et al.2009). Se- 05/04/2008 109 17-31 167 darknight vere quality cuts based on event rate after night sky IV 13/04/2008 97 17-22 147 darknight background suppression were applied to the sample; 29/04/2008 44 27-36 151 darknight 28.7hoursofgoodqualitydataremainedafterthesecuts, 03/05/2008 24 26-31 146 darknight out of which 11.5 hours were taken in 2007,8.7 hours in 04/05/2008 40 28-36 155 darknight 2008and 8.5 hours in 2009. More details about the final 05/05/2008 38 26-33 150 darknight 07/05/2008 40 28-36 153 darknight data set can be found in Table 1. For the signal study, 16/04/2009 93 17-27 133 2.8-3.7 a cut in the parameter size removed events with a to- 17/04/2009 103 17-28 151 1.6-2.4 talchargeless than 80 photo-electrons(phe) in the 2007 V 18/04/2009 126 17-28 168 0.7-1.7 data set, and 200phe in 2008and 2009data sets. In the 20/04/2009 73 19-35 171 0.9-1.4 latter case, this cut reduces the effect of the moon light. 21/04/2009 57 23-34 177 0.8-1.0 Finally, for the spectrum determination an additional 15/06/2009 57 24-35 125 0.8-2.1 cut in PMT DC current, namely above 2.5 µA, was ap- plied to the 2009 sample, in order to reduce systematics aPG 1553+113 data set from 2007 to 2009 used in this study. due to the moon light (Britzger et al. 2009), resulting Fromlefttoright: MAGICCycleofobservation,firstcolumn,and corresponding dates in dd/mm/yy, second column; effective time in 6.9hours of good quality data. For the conclusive ofobservationinminutes andzenithanglerangeindegrees, third steps of the analysis, Monte Carlo (MC) simulations of andfourthcolumn. Inthelasttwocolumns,therateoftheevents γ-like events were used. Hadronic background suppres- after the image cleaning, in Hz, and the mean DC current in the camera,inunitofµA,areshown. Thenightisconsideredasdark sionwasachievedusingtheRandomForest(RF)method night,iftheDCcurrentwhileobservinganextragalacticobjectis (Albert et al. 2008c), in which each event is assigned an lessthanindicatively1.2µA. additional parameter, the hadronness, which is related to the probability that the event is not γ-like. The RF The total Field of View (FoV) of the MAGIC camera method was also used in the energy estimation. The is 3.5◦, and the effective collection area is of the order threshold of the analysis was estimated to be 80GeV in of 105m2 at 200GeV for a source close to zenith. The 2007,150GeVin2008and160GeVin2009,asshownin incident light pulses are converted into analog signals, Table 2. transmittedviaopticalfibersanddigitisedby2GHzfast Duetochangesinthetelescopeperformance,thesigma analog to digital converters (FADCs). of the optical point-spread function (PSF) of 2007 and PG 1553+113 was observed with the MAGIC tele- 2008 was measured to be 13.0mm, while in 2009 it was scopefornearly19hoursin2005and2006(Albert et al. 14.9mm. To take all these differences into account, the 2007a); it was also the subject of a multiwavelength data were analysed separately, using dedicated sets of campaign carried out in July 2006 with optical, X-ray simulated data. and TeV γ ray telescopes (Albert et al. 2009). Here, we present the results of follow-up observations, performed 3. VHEγ RAYRESULTS for 14 hours in March-April 2007, for nearly 26 hours The 28.7 hours of good quality observations of in March-May 2008, some of those simultaneously with PG 1553+113 carried out between 2007 and 2009 re- other instruments (Aleksi´c et al. 2010), and for about sulted in a signal of 8.8σ of significance according to 24hoursinMarch-July2009,whichwerepartlytakenin eq.17 of Li & Ma (1983), obtained by combining the re- moderate light conditions (moon light). Unfortunately, sults from each year, listed in Table 2. The signal was both 2008 and 2009 observations were severely affected extractedbyanalyzingthe distributionofthe parameter 4 Aleksi´c et al. TABLE 3 PG 1553+113MEASUREDSPECTRA Year F >150GeV F >150GeV f0 Γ [cm−2 s−1] [Crab%] [cm−2 s−1 TeV−1] 2007 (1.40±0.38)·10−11 4% (1.1±0.3)·10−10 4.1±0.3 2008 (3.70±0.47)·10−11 11% (2.6±0.3)·10−10 4.3±0.4 2009 (1.63±0.45)·10−11 5% (1.3±0.2)·10−10 3.6±0.5 Note. — Spectra of the individual years of observations of PG 1553+113. From left to right: Year of MAGIC observations; Effective timeinhours;Integralfluxabove150GeVinunitsofcm−2 s−1 andCrabNebula%,normalizationfactorf0 inunitsofcm−2 s−1 TeV−1; inthelastcolumn,theΓindexobtainedbyfittingtheobserveddifferentialspectrumwithapowerlaw. Theerrorsreportedarestatistical only. Thesystematicuncertainty isestimatedtobe35%inthefluxleveland0.2inthepowerindex. alpha, related to the incoming direction of the primary for the 2007, 2008, and 2009 data sets respectively, ac- cosmic ray inducing the atmospheric shower. More de- cording to the significance of the signal. The 2008 data tails on the signal extraction with the alpha technique are consistent with the hypothesis of constant flux with can be found in Albert et al. 2008b. For the signal de- a probability of 50% (Figure 2, lower panel). tection, no cut in energy was applied. In2007,theobservedtimeseriesisconsistentwithcon- The significance of the signal was 5.8σ in 2007, 8.1σ stant flux (93% of probability). Nothing similar can be in 2008 and 4.2σ in 2009. Due to a large difference in concludedforthe2009observationssincethesignificance the energy thresholds and changes in the experimental of the signal is too low. conditions, the obtained fluxes cannot be compared di- In general, the high energy threshold of the analysis rectly. A detailed spectral analysis is necessary in order together with the weakness of the PG 1553+113 signal to study the source emission. and its very steep spectrum make any variability study at short timescale difficult and might have hidden the 3.1. Integral Flux detectionofanincreasedactivityonveryshorttimescale. In order to explore the VHE γ ray emission of 3.2. Differential Flux PG 1553+113 from each year, we compared the inte- gralflux above150GeV.Thisvalueisasafecompromise The differential spectra observed from PG 1553+113 taking into account the different energy thresholds. The by MAGIC each year from 2007 to 2009 are shown in final samples and the results of the spectral analysesare the left plot of Figure 3. shown in Table 3. As for other blazars, each spectrum can be well fitted The integral fluxes measured above 150GeV lie in the with a power law function of the form range of 4% to 11% of the Crab Nebula flux measured −Γ dF E by MAGIC (Albert et al. 2008a): the highest flux level =f ∗ (1) 0 is recorded in 2008 (0.11 Crab units1), a factor between dE (cid:18)200GeV(cid:19) twotothreelargercomparedtotheonemeasuredin2007 where f is the flux at 200GeV and Γ is the power law (0.04Crabunits)and2009(0.05Crabunits). Aconstant 0 index. The resulting indices are listed in the last col- fit to the data has a probability smaller than 0.001%. umnofTable3. Thesystematicuncertaintyisestimated Such changes in the flux level observedin PG 1553+113 to be 35% in the flux level and 0.2 in the power index are quite moderate in comparison to other monitored (Albert et al. 2008a), and is the sum of many contribu- TeV blazars. For example, in Mkn 421 a flux variations tions, mainly related to the use of MC simulations in- exceeding one order of magnitude have been observed stead of test beams. Thanks to the low energy thresh- (e.g. Fossati et al 2008). old of the analysisof 2007data, the correspondingspec- Adetailedstudyaboutpossiblefluxlevelvariationson trum has a measured point below 100GeV. In particu- short timescale was carried out with the limiting condi- lar, the measured differential energy flux at 98GeV is tion that the signal is not strong enough to allow for a (2.7 ± 0.3)·10−9cm−2s−1TeV−1, in agreement, within detailedsamplingonsub-daytimescale. Theupperpanel the errors, with the low energy point measured in 2005 ofFigure1displaysthelightcurveofPG1553+113mea- and 2006, (4.1 ± 1.2)·10−9cm−2s−1TeV−1 at 97 GeV. suredfrom 2007to 2009by MAGIC with a variable bin- The2008differentialenergyspectrummeasuredabove ning. For comparison, the daily flux levels measured in 150GeV has a slope of 4.3±0.4, while the slope of 2005and2006areshown,as extrapolatedfromthe pub- the spectrum determined with a partially simultaneous lished data (Albert et al. 2007a), and rescaled accord- sample taken during a multiwavelength campaign with ing to the power laws that interpolate the differential other instruments is 3.4±0.1 between 70 to 350GeV fluxes. Furthermore, the 2006 mean integral flux above (Aleksi´c et al. 2010). The different energy range char- 150GeVtakenduringthemultiwavelengthcampaignand acterizing the measurements fully accounts for this ap- reported in Albert et al. (2009) is shown. The former parentdisagreement: thespectralpointsmeasuredinthe data have not been used for the integral flux study, due range 150–350GeV are, in fact, in very good agreement. to very large uncertainties related to the extrapolation Finally, the 2009 differential spectrum is barely deter- procedure. We set 2-days, daily and monthly binning mined due to the limited signal. Except for the latter sample,whose significanceis ratherlow andcorrespond- 1TheCrabunitusedinthisworkisanarbitraryunitobtainedby ingerrorsnoticeablylarge,thepowerlawindicesdescrib- dividingtheintegralenergyfluxmeasuredaboveacertainthresh- old by the Crab Nebula flux measured above the same threshold ing the spectra are compatible. This indicates that the byMAGIC(Albertetal.2008a). shape of the emitted spectrum does not change, even if PG 1553+113: five years of observations with MAGIC 5 2005 2006 2007 2008 2009 -110 1 10 0 GeV)x -2] sm-1 5 VHE γ rays (MAGIC) 15 [c > F (E 0 53500 54000 54500 55000 d) 15 Optical (KVA+Tuorla) an y] R b mJ 10 F ( [ 5 -110 20 53500 54000 54500 55000 1 0.3 - 10 keV)x -2] s[erg cm-1 1105 X rays (Swift/XRT) F ( 5 53500 54000 54500 55000 -80 3 HE γ rays (Fermi/LAT) 1 1 GeV)x -2 ] sm-1 2 > c E [ F ( 1 0 53500 54000 54500 55000 55553535355333533543534505355535063565350073505735008455084500945059450004505045001450514500245052450034505345004450544500545055450064505645007450574500855058550095059500050505001051002052000500050000 Time [MJD] Fig.1.—MultiwavelengthlightcurveofPG1553+113 from2005tolate2009. Fromuppertolowerpanel: VHEγ raysabove150GeV measuredbyMAGIC(filledcircles),opticalfluxintheR-band(triangles),Xrays0.3–10keVflux(squares)andγ raysabove1GeV(open circles). Theerrorbarsreportedhave1σsignificance. Downwardtriangles,inlowerpanel,referto2sigmaupperlimitsonthesourceflux above1GeV. correct for the effects due to the finite energy resolution 17 (procedure called unfolding, Albert et al. (2007b)). The F (R band)[mJy]1165 ggscoeasotldes,atthghareetemdmeeesanpntitflaeumtxhoenemg(sittmhteeadslle)bmyvaetrahinaisbdisleoittueyrrmcseeieninsatsoitonanbysleesa.urglAy- 14 powerlawfitgivesthe valuesΓ =4.27±0.14forthe in- dex and f = (1.61 ± 0.14) · 10−10 s−1 cm−2 TeV−1 for 54540 54550 54560 54570 54580 54590 0 -1110 6 rthesepnoonrdminaglizpartoibonabfialicttyoro,fw7i7th%a. χT2h/edoinft=eg5ra.7l/fl9uaxndabcoovre- >150 GeV)x-2-1] s[cm 24 1so5Vr0bHGedEeVinpishtoahttoetnhisnetflererovamecltioocfno8sm%woitolhofgtithcheaelCEdraBisbtLaN.necTbeauskliaanrgfleuixna.tbo- F (E 0 account the EBL absorption assuming the background 54540 54550 54560 54570 54580 54590 model proposed in Dominguez et al. (2011) and a red- Time [MJD] shift z = 0.40, the intrinsic spectrum is compatible with Fig.2.— Zoom of Figure 1 around 2008 MAGIC observations a power law of index 3.09 ± 0.20, as drawn in Figure 3. (lower panel). In the upper panel the corresponding optical flux If we assume z = 0.45, the corresponding spectrum is hasbeensuperimposedforcomparison. compatible with a power law of index 2.91 ± 0.21. thetotalfluxshowshintsof(smallamplitude)variability, 4. PG1553+113ASSEENATOTHERWAVELENGTHS as also noted for other BL Lacs such as 1ES 1218+304 Figure1displaysthelightcurveofPG1553+113indif- (Acciari et al. 2010). ferent wavelengths. The MAGIC data shown cover five The right plot of Figure 3 shows the combined dif- cycles of TeV observations at energies above 150GeV. ferential spectrum of PG 1553+113 from 2007 to 2009, The time bins used are variable,as described in the pre- superimposed to the 2005-2006 spectrum measured by vious section. MAGIC(Γ=4.21±0.25,Albert et al.2007a). Thegray The simultaneous optical R-band data are outlined in band represents the systematic effect on the combined the second panel. These data are collected on a nightly spectrum resulting from the use of different methods to basisby the TuorlaObservatoryBlazarMonitoring Pro- 6 Aleksi´c et al. -1] TeV10-8 -1]TeV10-8 -2-1 s [cmdN/dE1100-1-90 -2-1 s[cmdN/dE1100-1-90 10-11 10-11 10-12 10-12 2007 10-13 2008 2009 PG 1553+113 2005-2009 mean spectrum 10-13 Crab Flux (Albert et al. 2008a) 10-14 PG 1553+113 2005-2009 mean spectrum deabsorbed Crab Flux (Albert et al. 2008a) 102 102 Energy [GeV] Energy [GeV] Fig.3.—Differential energy spectra fromPG1553+113. LeftFigure: comparisonbetween 2007, 2008 and2009 spectra. Right Figure: superimpositionof2005-2006spectrum,fromAlbertetal.(2007a),to2007-2009meanspectrumandcorrespondingdeabsorptionforz=0.4 usingtheEBLmodelofDominguezetal.(2011). Inbothfigures,thefitoftheCrabNebulaspectrummeasuredbyMAGIC(Albertetal. 2008a)issuperimposedforcomparison. gram2 (Takalo et al. 2007) using the KVA 35cm tele- Fermi data presented in this paper are restricted to the scope at La Palma and the Tuorla 1meter telescope in 1GeV-100GeV energy range and were collected from Finland. MJD 54682 (2008 August 4) to MJD 55200 (2010 Jan- For the third panel in Fig 1 we used the 14 uary 4) in survey mode. Swift pointed observations (Gehrels et al. 2004) of An unbinned analysis was performed to produce the PG1553+113performedfrom2005-04-20to 2010-02-05. light curve with the standard analysis tool gtlike, in- We summed the data collected on 5 and 7 July 2009 in cluded in the Science Tools software package (version orderto haveenoughstatisticsto obtainagoodspectral v09r21p00). P6 V11 DIFFUSE Instrument Response fit. The XRT data were processed with standard proce- Functions (IRFs) were used, which are a refinement to dures (xrtpipeline v0.12.6), filtering, and screening previous analyses reflecting improved understanding of criteria by using the Heasoft package (v6.11). We con- the point spreadfunction and effective area(Abdo et al. sider the data collected in photon counting mode, and 2011, in preparation). For this analysis, only photons thus only XRT event grades 0–12 were selected. belonging to the Diffuse class and located in a circular ◦ Source events were extracted from a circular region Region Of Interest (ROI) of 10 radius, centered at the with a radius of 20 pixels (1 pixel ∼ 2.36”), while back- positionof PG 1553+113,were selected. In addition, we ◦ groundeventswereextractedfromacircularregionwith excluded photons arriving from zenith angles > 105 to radiusof50pixelsawayfromthesourceregion. Someob- limitcontaminationfromEarthlimbγ rays,andphotons ◦ servations showed an averagecount rate of > 0.5 counts with rocking angle > 52 to avoid time intervals during s−1,thuspile-upcorrectionwasrequired. Inthatcasewe which Earth entered the LAT FoV. extracted the source events from an annular region with A separate analysis of the high energy emission in an inner radius from 2 to 7 pixels (depending on the each time bin was performed. All point sources in the source count rate and estimated by means of the PSF 1FGL within 15 degrees of PG 1553+113, including the fitting technique, see Moretti et al. 2005) and an outer source of interest itself, were considered in the analysis. radius of 30 pixels. We extracted background events Those within the 10 degree radius ROI were fitted with within an annular region centered on the source with a power law with spectral indices frozen to the values radii 70 and 120 pixels. Ancillary response files were obtainedfromthe likelihoodanalysisofthefulldataset, generated with xrtmkarf, and account for different ex- whilethosebeyond10degreeradiusROIhadtheirvalues traction regions, vignetting and PSF corrections. We frozen to those found in 1FGL. usedthelastspectralredistributionmatricesintheCali- Upper limits at 2σ confidence level (downward trian- bration database maintained by HEASARC. All spectra gles in Figures 1 and 4) were computed for time bins were rebinned with a minimum of 20 counts per energy with Test Statistics (TS)3 < 4, and were handled as in bin to allow χ2 fitting within XSPEC (v12.7.0; Arnaud the first Fermi/LAT catalog paper. The estimated sys- 1996). We fit the spectrum with an absorbed (model tematic uncertainty on the flux is 10% at 100MeV, 5% tbabs in Xspec) log parabola law (see e.g. Tramacere et at 500MeV, and 20% at 10GeV. al. 2007), with a neutral hydrogen column fixed to its As already noted, during the five years of monitoring Galactic value (3.65×1020 cm−2; Kalberla et al. 2005). the source generally showed a marginal activity in the The observed 0.3-10 keV fluxes obtained by Swift/XRT VHE γ ray band. The same behavior is followed by the in the different observations are reported in the third optical flux, whose variations are limited within a fac- panel of Fig. 1. tor of four, with a maximum flux reached in 2008 and In the lower panel, the Fermi/LAT light curve of a minimum value in 2009. A low emission in 2009 is PG 1553+113, computed in 10-day bins, is displayed. also registered at all the other wavelengths, Figure 4, 3 TS is 2 times the difference of the log(likelihood) with and 2 Moreinformationathttp://users.utu.fi/kani/1m/ withoutthesource(Mattoxetal.1996). PG 1553+113: five years of observations with MAGIC 7 -1110 3 0 GeV)x -2] sm-1 2 VHE γ rays (MAGIC) 15 [c 1 > E 0 F ( 54800 54850 54900 54950 55000 55050 55100 55150 Optical (KVA+Tuorla) d) 10 an y] b J R m F ( [ 5 -110 54800 54850 54900 54950 55000 55050 55100 55150 V)x1 ]s-1 6 X rays (Swift/XRT) F (0.3 - 10 ke -2 [erg cm 24 54800 54850 54900 54950 55000 55050 55100 55150 3 -810 HE γ rays (Fermi/LAT) eV)x ]s-1 2 1 G -2 m > c E [ 1 F ( 54800 54850 54900 54950 55000 55050 55100 55150 5555555545454545544554458548458458454858454858450485084518451845824528458345384508445548450584555845609455694570945579458094558945099455994500945509450194551945029555295503955539550495554955059555595506055560550705557055080555805509055590550005550055010555105502055520550305553055040555405505055550550615556155075515755108155815091559150015501501515151021521031531041541510515065607570858095905050505 Time [MJD] Fig.4.—ZoomofFigure1aroundtheCycleVobservationscarriedoutin2009. Fromuppertolowerpanel: VHEγraysabove150GeV measured byMAGIC (filled circles), optical flux inthe R-band (triangles), X rays 0.3–10keV flux (squares) and softγ rays above 1GeV (opencircles). Downwardtriangles,inlowerpanel,referto2sigmaupperlimitsonthesourcefluxabove1GeV. suggestingthatthe sourceenteredinalowactivitystate duringthatyear,withaminimumreachedfewdaysafter -2-1]sm 8 MtwFAeeiGgnuICroepot5ibcsasehlrovawantdsiotTnhsee.Vressimulutltoafnaeocuosrroeblsaetriovnatisotnusd.yTbhee- -11 [cGeV)x10 567 VthHeEopγtircaaylflfluuxx.abInovoer1d5e0rGtoeVincisrepalsoettsetdataisstaicfsu,nwcteiounseodf > 150 4 E for 2007/8/9 samples the daily light curve values; how- F ( 3 ever,sincetheopticalmeasurementshaveadifferenttime 2 coverage, in some cases we derived the mean VHE flux fromtwoormore consecutivedays. 2005and2006data, 1 fromAlbert et al.(2007a),wererejectedfromthisstudy, 0 due to the large uncertainty on the extrapolated flux in -1 the VHE band. The mean flux value from 2006 multi- 8 10 12 14 16 wavelengthcampaign,reportedinAlbert et al.(2009),is F [mJy] included. A linear relation among the two components Fig.5.— Correlation study between PG 1553+113 optical R– has a 74% probability, which suggests a correlation be- band flux and VHE γ ray integral flux above 150GeV observed tween these two extreme energy bands. This result is from2006to2009. in good agreement with the SSC model, which predicts larger than that observed in the TeV, optical and GeV a correlation between the synchrotron and the IC emis- bands. The different variability displayed by the syn- sion, related to the same electron population. Due to chrotron(X-ray)and inverse Compton (GeV-TeV) com- the poor simultaneity ofVHE data with the other wave- ponents seems to be somewhat in contrastwith the typ- lengths,thesamestudyhasnotbeenperformedinX-rays icalbehavior observedin TeV BL Lacs,showing,in gen- and soft γ rays. eral, a coordinated variability (e.g. Fossati et al 2008)4. The X-ray light curve shows a pronounced variability, However, the sparse sampling of the observations and in contrast to optical and very high energy bands. The X-ray flux spans an interval of about one order of mag- 4 Rare exceptions to this rule are the so called “orphan” TeV nitude (with maximum in 2005 and minimum in 2009), flares,e.g. Krawczynskietal.(2004) 8 Aleksi´c et al. thelackofatrulysimultaneousmonitoringpreventsany We modelthe SEDwiththe one-zoneSSC modelfully strong conclusion. In particular, no optical nor gamma- describedinMaraschi & Tavecchio(2003). Theemission raydataareavailableduringtheperiodofthemaximum zoneissupposedtobesphericalwithradiusR,inmotion X-ray flux in 2005. Coordinated multifrequency moni- with bulk Lorentz factor Γ at an angle θ with respect to toring is necessary to further investigate this important the line ofsight. Special relativistic effects are described issue. by the relativistic Doppler factor, δ =[Γ(1−βcosθ)]−1. In a dedicated paper (Abdo et al. 2010b), the The energy distribution of the relativistic emitting elec- Fermi/LAT collaboration reported the analysis of the tronsisdescribedbyasmoothedbrokenpowerlawfunc- first year of PG 1553+113 data. They found that dur- tion, with limits γ and γ and break at γ . To cal- min max b ing the monitoring period the emission above 200MeV culate the SSC emission, we use the full Klein–Nishina is almost steady. This is in contrast with the be- cross section. havior of the source at higher energies (> 1GeV). Given the large variations of the X-ray synchrotron In our analysis of 2008 and 2009 LAT data (drawn flux, we decided to use the average level of the syn- in lower panel of Figure 1), in fact, a steady emis- chrotronbumpasmeasuredbyXRT,includingalsoASM sion above 1GeV has a probability smaller than 0.1% andBATfluxestoconstrainthemodel. Thecorrespond- and is ruled out. The lowest flux, observed in April inginputparametersarelistedinTable4. Wealsoreport 2009, has a value (0.5 ± 0.3) · 10−8 cm−2s−1, while the derived powers carried by the different components, the highest flux, detected in August 2008, has a value relativistic electrons, P , magnetic field, P , and pro- e B (2.8 ± 0.6) · 10−8 cm−2s−1, more than 5 times higher. tons, P , (assuming a composition of one cold proton p In our case, we are looking only at the upper edge of perrelativisticelectron)andthetotalradiativeluminos- LATband,probablyclosetotheICpeak,whiletheinte- ity L ≃L /δ2. r obs gralfluxabove200MeVreportedbyAbdo et al.(2010b) In order to investigate the role of different parameters is dominated by lower energies, due to larger statis- inthe model, we haveexploredtheir variationas a func- tics. Therefore, we conclude that while at low energies tion of the intensity of the synchrotron peak. To do so, (200MeV-1GeV)theICcontinuumshowsonlymarginal we have modeled the SED considering the two extreme variability,this is notthe caseinthe vicinity ofthe peak states of the synchrotron peak described above, respec- (some GeVs). In fact, small variability at GeV energies tively. For the SSC peak, instead, we have fixed the is a common feature of HBLs (Abdo et al. 2010a). VHE data. The two curves representing the models are superimposed in light gray to Figure 6. The parameters 5. MODELINGTHESED obtained, listed in the last two columns of Table 4, are In Figure 6, we assembled the SED of PG 1553+113 quite similar to the ones obtained when considering the using historical data and the MAGIC spectra described averagelevelofthesynchrotronbump,exceptforthetwo above. Open black squares displaying radio-optical variables B and K, the magnetic field and the electron data are from NED5. In the optical band, we also density. Indeed, these two parameters regulate the rel- show (red diamonds) the KVA minimum and maximum ative importance of synchrotron and SSC components. flux measured in the period covered by MAGIC 2005- The state characterized by a low synchrotron emission 2009 observations together with optical-UV fluxes from has larger B and smaller K values with respect to the Swift/UVOT(filledblacktriangles,fromTavecchio et al. mean state modeled above. Conversely, the high syn- 2010). FortheX-raydata,twoSwift/XRTspectrataken chrotron emission state has smaller values of B and a in 2005 (high flux state, red crosses, and intermediate larger K. state, black asterisks, from Tavecchio et al. 2010) are Finally,acomparisonwiththeSEDmodelobtainedin given, and a Suzaku spectrum taken in 2006 (continu- the multiwavelength campaign reported in Aleksi´c et al. ous red line, from Reimer et al. 2008). In addition, the (2010),revealsthatthe parametersusedforbuilding the average 15–150 keV flux measured by Swift/BAT dur- two models are quite similar. The major differences are ingthe first54monthsofsurvey(Cusumano et al.2010) thevalueoftheDopplerfactor,whichinourmodelisrel- is shown (black star), and the averageRXTE/ASM flux atively higher (δ=35) thanin the previousone (δ=23), betweenMarch1andMay31,2008(smallblacksquare), and that of the magnetic field (0.5G instead of 0.7G). from quick–look results provided by the RXTE/ASM ThisdifferenceismainlyduetothehigherSSCpeakfre- team6. quency that we find in our data, better defined by the The green triangles correspond to the LAT spectrum combined LAT and MAGIC spectra. averaged over ∼200 days (2008 August-2009 February) The derived value of the total jet power, P = jet from Abdo et al. (2010b). As discussed in that paper, P + P + P = 4 × 1044 erg/s, is consistent with e B p the flux above 200 MeV is rather stable, showing very the typical values inferred modelling similar sources small variability over the entire period of LAT observa- (e.g. Ghisellini et al. 2011). We use a relatively large tions. It is likely that the variability observed at the minimum electron Lorentz factor γ ∼ 103 in or- min highest energies is not important in determining the av- der to reproduce the hard MeV-GeV continuum tracked eraged spectrum due to the limited statistics. by LAT (photon index Γ = 1.68 ± 0.03). The For MAGIC, we report the 2005-2006 and 2007-2009 high value of γ implies that, as commonly derived min observedspectra(filledcircles)andthesamespectracor- in TeV BL Lacs, the relativistic electrons (and the mag- rectedfor the absorptionby the EBL using the modelof netic field, almost in equipartition) carry more power Dominguez et al. (2011) (red open circles). than the cold proton component. Another characteris- tic that PG 1553+113 shares with the other TeV BL 5 http://nedwww.ipac.caltech.edu/ Lacs is that the total luminosity L is larger than the 6 http://xte.mit.edu/asmlc/ r PG 1553+113: five years of observations with MAGIC 9 variability within a factor 4. A clear variability is seen TABLE 4 Inputmodelparametersforthe modelsshown in Fig.6 in the X-rays and γ rays above 1GeV. The latter out- come,exploringtheenergiesclosetotheICpeak,isonly Parameter Valuemean Valuemax Valuemin apparently in contradiction with previous results stat- γmin [103] 2.5 1 5 ing a quite stable spectrum for this source in the soft γb [104] 3.2 3 1.3 (> 200 MeV) γ–ray band (Abdo et al. 2010b). The dif- γmax [105] 2.2 5.2 4.1 ferent energy thresholds used in the two studies can, in n1 2.0 2.0 2.0 fact, explain very well the discrepancy, as discussed in n2 4.0 3.75 3.55 the paper. B [G] 0.5 0.8 0.2 K [103 cm−3] 5.35 3.8 25 Finally, for the study of the spectral energy distri- R [1016 cm] 1 1 1 bution, the mean differential spectrum measured by δ 35 35 35 MAGICwascombinedwithhistoricaldataatotherwave- Pe [1044 erg/s] 2.2 lengths. Due to the large variations observed in X-rays PB [1044 erg/s] 1.5 and characterizing the synchrotron peak, we decided to Pp [1044 erg/s] 0.34 use for the SED modeling the high energy bump, and Lr [1044 erg/s] 6.3 the average level of the low energy bump. A more pre- cise model requires coupling the VHE γ ray part of the aWe listforthethree differentmodelsplotted inFig.6themin- spectrumwith simultaneouscoverageofthe synchrotron imum,break andmaximum Lorentz factors and thelow and high peak,inparticularatoptical-X-rayenergies. Aninterest- energyslopeoftheelectronenergydistribution,themagneticfield ingfeatureofPG1553+113isthenarrownessoftheSSC intensity,theelectrondensity,theradiusoftheemittingregionand itsDopplerfactor. Fortheaveragemodelwealsogivethederived peak derived fromthe LAT and MAGIC spectra, imply- power carriedby electrons, magnetic field, protons (assumingone ingarelativelylargevalueoftheminimumLorentzfactor coldprotonperemittingrelativisticelectron) andthetotal radia- ofthe emitting electrons,2.5×103. This is alsorequired tiveluminosity. by other HBLs with hard GeV spectra (e.g. Tavecchio power supplied by electrons, magnetic field and protons. et al. 2010). As discussed in Celotti & Ghisellini (2008), this implies TheMAGICstereosystem,withitsincreasedsensitiv- that either only a small fraction of leptons is acceler- ity and low energy threshold, is the suitable instrument ated at relativistic energies (leaving a reservoir of cold to further investigate eventual daily scale TeV variabil- pairs and/or protons) or that the jet is dissipating a ity, as well as to provide a good differential spectrum large fraction of its power as radiation, eventually lead- determination below 100GeV. ing to the deceleration of the flow, as in fact observed at VLBI scales (Piner et al. 2010) and envisaged in the models of structured jets (Georganopoulos & Kazanas ACKNOWLEDGEMENTS 2003; Ghisellini et al. 2005). We would like to thank the Instituto de Astrof´ısica 6. CONCLUSIONS de Canarias for the excellent working conditions at the Observatorio del Roque de los Muchachos in La Palma. In this paper, we presented the analysis of three years The support of the German BMBF and MPG, the Ital- ofVHEγ raydataofPG1553+113collectedbyMAGIC ian INFN, the Swiss National Fund SNF, and the Span- from 2007 to 2009. The data set was divided into indi- ish MICINN is gratefully acknowledged. This work was vidual years, and a significant signal was found in every also supported by the Marie Curie program, by the sample, confirming PG 1553+113 as a stable presence CPAN CSD2007-00042 and MultiDark CSD2009-00064 in the VHE sky. The overall flux above 150GeV from projects of the Spanish Consolider-Ingenio 2010 pro- 2007 to 2009 shows only a modest variability on yearly gramme, by grant DO02-353 of the Bulgaria n NSF, by time-scale, within a factor 3, corresponding to a varia- grant 127740 of the Academy of Finland, by the YIP of tion between 4% to 11% of the Crab Nebula flux. No theHelmholtz Gemeinschaft,bytheDFGClusterofEx- clear variability on smaller time scales is evident in the cellence “Origin and Structure of the Universe”, by the sample. DFG Collaborative Research Centers SFB823/C4 and For the spectral analysis, the data set was combined SFB876/C3, and by the Polish MNiSzW grant 745/N- withpreviousobservationscarriedoutbyMAGICduring HESS-MAGIC/2010/0. the firsttwoCyclesofoperations,from2005to 2006,for The Fermi LAT Collaborationacknowledges generous a total of five years of monitoring. This sample was ex- ongoing support from a number of agencies and insti- cludedfromthetemporalstudyduetoverylargesystem- tutesthathavesupportedboththedevelopmentandthe aticsrelatedtothefluxextrapolationprocedure. Despite operation of the LAT as well as scientific data analysis. the hints of variability on the flux level, the differential These include the National Aeronautics and Space Ad- flux from each year is in very good agreement with a ministrationandtheDepartmentofEnergyintheUnited power law of constant index 4.27 ± 0.14. This behavior States, the Commissariat `a l’Energie Atomique and the has been already observed in other blazars, such as the Centre National de la Recherche Scientifique / Institut HBL 1ES 1218+304(Acciari et al. 2010). National de Physique Nucl´eaire et de Physique des Par- PG1553+113wasalsomonitoredinoptical,X-rayand ticules in France, the Agenzia Spaziale Italiana and the softγ rayfrequencies,butonlythe formerdatacouldbe Istituto Nazionale di Fisica Nucleare in Italy, the Min- used for correlation studies thanks to the large timing istryofEducation,Culture,Sports,ScienceandTechnol- coverage. Interestingly, a hint of correlation with prob- ogy(MEXT),HighEnergyAcceleratorResearchOrgani- ability of 74% was found between MAGIC and R-band zation(KEK)andJapanAerospaceExplorationAgency optical flux levels, which in turn shows only a modest (JAXA)inJapan,andtheK.A.WallenbergFoundation, 10 Aleksi´c et al. Fig. 6.— SED of PG 1553+113. Open black squares are radio-optical data from NED, red diamonds represent the KVA minimum and maximum flux measured in the period covered by MAGIC observations, together with optical-UV fluxes from Swift/UVOT (filled black triangles). For the X-ray data, two Swift/XRT spectra taken in 2005 (high flux state, red crosses, and intermediate state, black asterisks) are given, and a Suzaku spectrum taken in 2006 (continuous red line). The average 15–150 keV flux measured by Swift/BAT (black star) and the average RXTE/ASM flux between March 1 and May 31, 2008 (small black square), are also represented. The LAT spectrum averaged over ∼200 days (2008 August-2009 February) is given (green filled triangles). For MAGIC, we display the 2005-2006 and2007-2009 observedspectra(bluefilledcircles)andthesamespectracorrectedfortheabsorptionbytheEBL(redopencircles). The average SED is modeled with a one-zone SSC model (continuous black line). Alternative SED are superimposed in light gray. Detailed referencesareaddressedinthetext. the Swedish ResearchCouncil and the Swedish National tituto Nazionale di Astrofisica in Italy and the Centre Space Board in Sweden. National d’E´tudes Spatiales in France. Additional support for science analysis during the op- erations phase is gratefully acknowledged from the Is- REFERENCES Abdo,A.A.,etal.2010a, ApJ,722,520 Albert,J.,etal.2008c, Nucl.Instr.Meth.A588,424 Abdo,A.A.,etal.2010b, ApJ,708,1310 Albert,J.,etal.2007a, ApJLetters,654,L119 Abdo,A.A.,etal.2009, ApJ,700,597 Albert,J.,etal.,2007b,Nucl.Instr.Meth.A583,494 Acciari,V.A.,etal.2010, ApJ,709,L163 Aleksi´c,J.,etal.2010,A&A515,A76 Aharonian,F.,etal.2006a, A&A,448,L19 Aliu,E.,etal.2009,ApJ,30,293 Aharonian,F.,etal.2006b, Nature440,1018 Arnaud,K.A.1996,ASPC,101,17 Albert,J.,etal.2009,A&A493,467 Baixeras,C.,etal.2004,Nucl.Instrum.Meth.A518,188 Albert,J.,etal.2008a, ApJ,674,1037 Albert,J.,et.al.2008b,ApJ,685,L23