# Incremental Fermi Large Area Telescope Fourth Source Catalog

Soheila Abdollahi, Fabio Acero, Luca Baldini, Jean Ballet, Denis Bastieri, Ronaldo Bellazzini, Bijan Berenji, Alessandra Berretta, Elisabetta Bissaldi, Roger D. Blandford, et al.

► To cite this version:

Soheila Abdollahi, Fabio Acero, Luca Baldini, Jean Ballet, Denis Bastieri, et al. Incremental Fermi Large Area Telescope Fourth Source Catalog. Astrophys.J.Supp., 2022, 260 (2), pp.53. 10.3847/1538-4365/ac6751. hal-03589174

HAL Id: hal-03589174

https://hal.science/hal-03589174

Submitted on 11 Aug 2022

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
<figure><img src="image_2.png" /></figure>

Incremental Fermi Large Area Telescope Fourth Source Catalog

S. Abdollahi¹, F. Acero² ID, L. Baldini³ ID, J. Ballet² ID, D. Bastieri⁴,⁵,⁶ ID, R. Bellazzini⁷ ID, B. Berenji⁸, A. Berretta⁹, E. Bissaldi¹⁰,¹¹ ID, R. D. Blandford¹² ID, E. Bloom¹², R. Bonino¹³,¹⁴ ID, A. Brill¹⁵ ID, R. J. Britto¹⁶ ID, P. Bruel¹⁷ ID, T. H. Burnett¹⁸ ID, S. Buson¹⁹, R. A. Cameron¹² ID, R. Caputo²⁰, P. A. Caraveo²¹ ID, D. Castro²⁰,²² ID, S. Chaty² ID, C. C. Cheung²³ ID, G. Chiaro²¹, N. Cibrario¹³, S. Cipriņi²⁴,²⁵ ID, J. Coronado-Blázquez²⁶,²⁷, M. Crnogorcevic²⁸, S. Cutini²⁹ ID, F. D'Ammando³⁰ ID, S. De Gaetano¹⁰,¹¹ ID, S. W. Digel¹², N. Di Lalla¹² ID, F. Dirirsa³¹, L. Di Venere¹⁰,¹¹ ID, A. Domínguez³² ID, V. Fallah Ramazani³³ ID, S. J. Fegan¹⁷, E. C. Ferrara²⁰,²⁸,³⁴ ID, A. Fiori³⁵ ID, H. Fleischhack²⁰,³⁴,³⁶, A. Franckowiak³³, Y. Fukazawa³⁷ ID, S. Funk³⁸ ID, P. Fusco¹⁰,¹¹ ID, G. Galanti²¹, V. Gammaldi²⁶,²⁷ ID, F. Gargano¹¹ ID, S. Garrappa³⁹ ID, D. Gasparrini²⁴,²⁵ ID, F. Giacchino²⁴,²⁵, N. Giglietto¹⁰,¹¹, F. Giordano¹⁰,¹¹, M. Giroletti³⁰ ID, T. Glanzman¹² ID, D. Green⁴⁰ ID, I. A. Grenier², M.-H. Grondin⁴¹ ID, L. Guillemot⁴²,⁴³ ID, S. Guiriec²⁰,⁴⁴ ID, M. Gustafsson⁴⁵, A. K. Harding⁴⁶ ID, E. Hays²⁰ ID, J. W. Hewitt⁴⁷ ID, D. Horan¹⁷, X. Hou⁴⁸,⁴⁹,⁵⁰ ID, G. Jóhannesson⁵¹,⁵² ID, C. Karwin⁵³, T. Kayanoki³⁷, M. Kerr²³ ID, M. Kuss⁷ ID, D. Landriu², S. Larsson⁵⁴,⁵⁵,⁵⁶ ID, L. Latronico¹³ ID, M. Lemoine-Goumard⁴¹ ID, J. Li⁵⁷ ID, I. Liodakis⁵⁸, F. Longo⁵⁹,⁶⁰ ID, F. Loparco¹⁰,¹¹ ID, B. Lott⁴¹ ID, P. Lubrano²⁹ ID, S. Maldera¹³ ID, D. Malyshev³⁸, A. Manfreda³ ID, G. Martí-Devesa⁶¹, M. N. Mazziotta¹¹ ID, I. Mereu⁹,²⁹ ID, M. Meyer³⁸, P. F. Michelson¹², N. Mirabal²⁰,⁶² ID, W. Mitthumsiri⁶³ ID, T. Mizuno⁶⁴ ID, A. A. Moiseev²⁸,³⁴, M. E. Monzani¹² ID, A. Morselli²⁴ ID, I. V. Moskalenko¹² ID, M. Negro³⁴,⁶², E. Nuss⁶⁵ ID, N. Omodei¹² ID, M. Orienti³⁰, E. Orlando¹²,⁶⁶, D. Paneque⁴⁰, Z. Pei⁵, J. S. Perkins²⁰, M. Persic⁵⁹,⁶⁷ ID, M. Pesce-Rollins⁷ ID, V. Petrosian¹² ID, R. Pillera¹⁰,¹¹ ID, H. Poon³⁷ ID, T. A. Porter¹² ID, G. Principe³⁰,⁵⁹,⁶⁰ ID, S. Rainò¹⁰,¹¹ ID, R. Rando⁴,⁵,⁶ ID, B. Rani²⁰,⁶⁸,⁶⁹ ID, M. Razzano³ ID, S. Razzaque⁷⁰ ID, A. Reimer¹²,⁶¹ ID, O. Reimer⁶¹ ID, T. Reposeur⁴¹, M. Sánchez-Conde²⁶,²⁷ ID, P. M. Saz Parkinson⁷¹,⁷²,⁷³ ID, L. Scotton⁶⁵, D. Serini¹¹ ID, C. Sgrò⁷ ID, E. J. Siskind⁷⁴, D. A. Smith⁴¹,⁷⁵ ID, G. Spandre⁷, P. Spinelli¹⁰,¹¹, K. Sueoka³⁷, D. J. Suson⁷⁶ ID, H. Tajima¹²,⁷⁷, D. Tak²⁰,⁷⁸ ID, J. B. Thayer¹², D. J. Thompson²⁰ ID, D. F. Torres⁷⁹,⁸⁰,⁸¹ ID, E. Troja²⁰,²⁸ ID, J. Valverde²⁰,⁶² ID, K. Wood⁸²,⁸⁴, and G. Zaharijas⁸³

¹ IRAP, Université de Toulouse, CNRS, UPS, CNES, F-31028 Toulouse, France

² Université Paris Saclay and Université Paris Cité, CEA, CNRS, AIM, F-91191 Gif-sur-Yvette, France; jean.ballet@cea.fr

³ Università di Pisa and Istituto Nazionale di Fisica Nucleare, Sezione di Pisa I-56127 Pisa, Italy

⁴ Istituto Nazionale di Fisica Nucleare, Sezione di Padova, I-35131 Padova, Italy

⁵ Dipartimento di Fisica e Astronomia “G. Galilei”, Università di Padova, I-35131 Padova, Italy

⁶ Center for Space Studies and Activities “G. Colombo”, University of Padova, Via Venezia 15, I-35131 Padova, Italy

⁷ Istituto Nazionale di Fisica Nucleare, Sezione di Pisa, I-56127 Pisa, Italy

⁸ California State University, Los Angeles, Department of Physics and Astronomy, Los Angeles, CA 90032, USA

⁹ Dipartimento di Fisica, Università degli Studi di Perugia, I-06123 Perugia, Italy

¹⁰ Dipartimento di Fisica “M. Merlin” dell’Università e del Politecnico di Bari, via Amendola 173, I-70126 Bari, Italy

¹¹ Istituto Nazionale di Fisica Nucleare, Sezione di Bari, I-70126 Bari, Italy

¹² W.W. Hansen Experimental Physics Laboratory, Kavli Institute for Particle Astrophysics and Cosmology, Department of Physics and SLAC National Accelerator Laboratory, Stanford University, Stanford, CA 94305, USA

¹³ Istituto Nazionale di Fisica Nucleare, Sezione di Torino, I-10125 Torino, Italy

¹⁴ Dipartimento di Fisica, Università degli Studi di Torino, I-10125 Torino, Italy

¹⁵ NASA Postdoctoral Program Fellow, NASA Goddard Space Flight Center, Greenbelt, MD 20771, USA

¹⁶ Department of Physics, University of the Free State, P.O. Box 339, Bloemfontein 9300, South Africa

¹⁷ Laboratoire Leprince-Ringuet, École polytechnique, CNRS/IN2P3, F-91128 Palaiseau, France; Philippe.Bruel@llr.in2p3.fr

¹⁸ Department of Physics, University of Washington, Seattle, WA 98195-1560, USA; tburnett@u.washington.edu

¹⁹ Institut für Theoretische Physik and Astrophysik, Universität Würzburg, D-97074 Würzburg, Germany

²⁰ NASA Goddard Space Flight Center, Greenbelt, MD 20771, USA

²¹ INAF-Istituto di Astrofisica Spaziale e Fisica Cosmica Milano, via E. Bassini 15, I-20133 Milano, Italy

²² Harvard-Smithsonian Center for Astrophysics, Cambridge, MA 02138, USA

²³ Space Science Division, Naval Research Laboratory, Washington, DC 20375-5352, USA

²⁴ Istituto Nazionale di Fisica Nucleare, Sezione di Roma “Tor Vergata”, I-00133 Roma, Italy

²⁵ Space Science Data Center—Agenzia Spaziale Italiana, Via del Politecnico, snc, I-00133, Roma, Italy

²⁶ Instituto de Física Teórica UAM/CSIC, Universidad Autónoma de Madrid, E-28049 Madrid, Spain

²⁷ Departamento de Física Teórica, Universidad Autónoma de Madrid, E-28049 Madrid, Spain

²⁸ Department of Astronomy, University of Maryland, College Park, MD 20742, USA

²⁹ Istituto Nazionale di Fisica Nucleare, Sezione di Perugia, I-06123 Perugia, Italy

³⁰ INAF Istituto di Radioastronomia, I-40129 Bologna, Italy

³¹ Astronomy and Astrophysics Research Development Department, Entoto Observatory and Research Center, Ethiopian Space Science and Technology Institute, Ethiopia

³² Grupo de Altas Energías, Universidad Complutense de Madrid, E-28040 Madrid, Spain

³³ Ruhr University Bochum, Faculty of Physics and Astronomy, Astronomical Institute (AIRUB), D-44780 Bochum, Germany

³⁴ Center for Research and Exploration in Space Science and Technology (CRESST) and NASA Goddard Space Flight Center, Greenbelt, MD 20771, USA

³⁵ Dipartimento di Fisica “Enrico Fermi”, Università di Pisa, Pisa I-56127, Italy

³⁶ Catholic University of America, Washington, DC 20064, USA

³⁷ Department of Physical Sciences, Hiroshima University, Higashi-Hiroshima, Hiroshima 739-8526, Japan

³⁸ Friedrich-Alexander Universität Erlangen-Nürnberg, Erlangen Centre for Astroparticle Physics, Erwin-Rommel-Str. 1, 91058 Erlangen, Germany

³⁹ Deutsches Elektronen Synchrotron DESY, D-15738 Zeuthen, Germany

⁴⁰ Max-Planck-Institut für Physik, D-80805 München, Germany

⁴¹ Université Bordeaux, CNRS, LP2I Bordeaux, UMR 5797, F-33170 Gradignan, France; lott@cenbg.in2p3.fr

⁴² Laboratoire de Physique et Chimie de l'Environnement et de l'Espace—Université d'Orléans/CNRS, F-45071 Orléans Cedex 02, France
⁴³ Station de radioastronomie de Nançay, Observatoire de Paris, CNRS/INSU, F-18330 Nançay, France

⁴⁴ The George Washington University, Department of Physics, 725 21st Street, NW, Washington, DC 20052, USA

⁴⁵ Georg-August University Göttingen, Institute for Theoretical Physics—Faculty of Physics, Friedrich-Hund-Platz 1, D-37077 Göttingen, Germany

⁴⁶ Los Alamos National Laboratory, Los Alamos, NM 87545, USA

⁴⁷ University of North Florida, Department of Physics, 1 UNF Drive, Jacksonville, FL 32224, USA

⁴⁸ Yunnan Observatories, Chinese Academy of Sciences, 396 Yangfangwang, Guandu District, Kunming 650216, People's Republic of China

⁴⁹ Key Laboratory for the Structure and Evolution of Celestial Objects, Chinese Academy of Sciences, 396 Yangfangwang, Guandu District, Kunming 650216, People's Republic of China

⁵⁰ Center for Astronomical Mega-Science, Chinese Academy of Sciences, 20A Datun Road, Chaoyang District, Beijing 100012, People's Republic of China

⁵¹ Science Institute, University of Iceland, IS-107 Reykjavik, Iceland

⁵² Nordita, Royal Institute of Technology and Stockholm University, Roslagstullsbacken 23, SE-106 91 Stockholm, Sweden

⁵³ Department of Physics and Astronomy, Clemson University, Kinard Lab of Physics, Clemson, SC 29634-0978, USA

⁵⁴ Department of Physics, KTH Royal Institute of Technology, AlbaNova, SE-106 91 Stockholm, Sweden

⁵⁵ The Oskar Klein Centre for Cosmoparticle Physics, AlbaNova, SE-106 91 Stockholm, Sweden

⁵⁶ School of Education, Health and Social Studies, Natural Science, Dalarna University, SE-791 88 Falun, Sweden

⁵⁷ Department of Astronomy, School of Physical Sciences, University of Science and Technology of China, Hefei, Anhui 230026, People's Republic of China

⁵⁸ Finnish Centre for Astronomy with ESO (FINCA), University of Turku, FI-21500 Piikiiö, Finland

⁵⁹ Istituto Nazionale di Fisica Nucleare, Sezione di Trieste, I-34127 Trieste, Italy

⁶⁰ Dipartimento di Fisica, Università di Trieste, I-34127 Trieste, Italy

⁶¹ Institut für Astro- und Teilchenphysik, Leopold-Franzens-Universität Innsbruck, A-6020 Innsbruck, Austria

⁶² Department of Physics and Center for Space Sciences and Technology, University of Maryland Baltimore County, Baltimore, MD 21250, USA

⁶³ Department of Physics, Faculty of Science, Mahidol University, Bangkok 10400, Thailand

⁶⁴ Hiroshima Astrophysical Science Center, Hiroshima University, Higashi-Hiroshima, Hiroshima 739-8526, Japan

⁶⁵ Laboratoire Univers et Particules de Montpellier, Université Montpellier, CNRS/IN2P3, F-34095 Montpellier, France

⁶⁶ Istituto Nazionale di Fisica Nucleare, Sezione di Trieste, Università di Trieste, I-34127 Trieste, Italy

⁶⁷ Osservatorio Astronomico di Trieste, Istituto Nazionale di Astrofisica, I-34143 Trieste, Italy

⁶⁸ Korea Astronomy and Space Science Institute, 776 Daedeokdae-ro, Yuseong-gu, Daejeon 30455, Republic of Korea

⁶⁹ Department of Physics, American University, Washington, DC 20016, USA

⁷⁰ Centre for Astro-Particle Physics (CAPP) and Department of Physics, University of Johannesburg, PO Box 524, Auckland Park 2006, South Africa

⁷¹ Santa Cruz Institute for Particle Physics, Department of Physics and Department of Astronomy and Astrophysics, University of California at Santa Cruz, Santa Cruz, CA 95064, USA

⁷² Department of Physics, The University of Hong Kong, Pokfulam Road, Hong Kong, People's Republic of China

⁷³ Laboratory for Space Research, The University of Hong Kong, Hong Kong, People's Republic of China

⁷⁴ NYCB Real-Time Computing Inc., Lattingtown, NY 11560-1025, USA

⁷⁵ Laboratoire d'Astrophysique de Bordeaux, Université de Bordeaux, CNRS, B18N, allée Geoffroy Saint-Hilaire, F-33615 Pessac, France

⁷⁶ Purdue University Northwest, Hammond, IN 46323, USA

⁷⁷ Solar-Terrestrial Environment Laboratory, Nagoya University, Nagoya 464-8601, Japan

⁷⁸ Department of Physics, University of Maryland, College Park, MD 20742, USA

⁷⁹ Institute of Space Sciences (ICE, CSIC), Campus UAB, Carrer de Magrans s/n, E-08193 Barcelona, Spain

⁸⁰ Institut d'Estudis Espacials de Catalunya (IEEC), E-08034 Barcelona, Spain

⁸¹ Institució Catalana de Recerca i Estudis Avançats (ICREA), E-08010 Barcelona, Spain

⁸² Praxis Inc., Alexandria, VA 22303, USA

⁸³ Center for Astrophysics and Cosmology, University of Nova Gorica, Nova Gorica, Slovenia

Received 2022 January 26; revised 2022 April 11; accepted 2022 April 12; published 2022 June 30

## Abstract

We present an incremental version (4FGL-DR3, for Data Release 3) of the fourth Fermi Large Area Telescope (LAT) catalog of γ-ray sources. Based on the first 12 years of science data in the energy range from 50 MeV to 1 TeV, it contains 6658 sources. The analysis improves on that used for the 4FGL catalog over eight years of data: more sources are fit with curved spectra, we introduce a more robust spectral parameterization for pulsars, and we extend the spectral points to 1 TeV. The spectral parameters, spectral energy distributions, and associations are updated for all sources. Light curves are rebuilt for all sources with 1 yr intervals (not 2 month intervals). Among the 5064 original 4FGL sources, 16 were deleted, 112 are formally below the detection threshold over 12 yr (but are kept in the list), while 74 are newly associated, 10 have an improved association, and seven associations were withdrawn. Pulsars are split explicitly between young and millisecond pulsars. Pulsars and binaries newly detected in LAT sources, as well as more than 100 newly classified blazars, are reported. We add three extended sources and 1607 new point sources, mostly just above the detection threshold, among which eight are considered identified, and 699 have a plausible counterpart at other wavelengths. We discuss the degree-scale residuals to the global sky model and clusters of soft unassociated point sources close to the Galactic plane, which are possibly related to limitations of the interstellar emission model and missing extended sources.

Unified Astronomy Thesaurus concepts: Gamma-ray sources (633); Gamma-ray astronomy (628); Blazars (164); Pulsars (1306); Sky surveys (1464); Celestial objects catalogs (212)

## 1. Introduction

The Fermi Large Area Telescope (LAT) has been surveying the high-energy γ-ray sky since 2008 (Atwood et al. 2009), and the LAT Collaboration has published a succession of source catalogs based on comprehensive analyses of LAT data. The

<figure><img src="image_3.png" /></figure>

<footnote n="84">⁸⁴ Resident at Naval Research Laboratory, Washington, DC 20375, USA.</footnote>

<footnote>Original content from this work may be used under the terms of the Creative Commons Attribution 4.0 licence. Any further distribution of this work must maintain attribution to the author(s) and the title of the work, journal citation and DOI.</footnote>
fourth source catalog (4FGL, where FGL stands for Fermi Gamma-ray LAT) was derived from the analysis of the first 8 yr of LAT science data and contained 5064 sources (Abdollahi et al. 2020, hereafter the 4FGL paper, or simply 4FGL). Every FGL catalog until and including 4FGL had used a new analysis method, new calibrations, a new diffuse model, and even a new reconstruction of the events themselves (like Pass 8 for 4FGL). Therefore each successive version of the catalog was largely independent of the previous one. Results for different versions were compared only after the fact.

Since the development of the 4FGL, the data have been stable, no new model for the Galactic interstellar emission was built, and the analysis method has not evolved much. So we decided to change our approach and start building incremental 4FGL versions every two years, until one major analysis or data component changes. The Data Release 2 (DR2; Ballet et al. 2020) covered 10 yr of data and contained 723 additional sources. We refer to the original 4FGL sources as Data Release 1 (DR1).

In this paper, we describe the incremental procedure and apply it to 12 yr of data, creating the Data Release 3 (4FGL-DR3 or DR3 for short). The 4FGL paper remains the reference for the detailed methodology. DR2 was strictly incremental (all DR1 sources are included and the analysis procedure was exactly the same). Meanwhile we have improved a few analysis steps, and recent publications have revealed a few new extended sources. Therefore we have relaxed this condition in DR3, to allow for deleting a few sources, and accept a few improved localizations.

All the information present in the DR1 catalog is updated, except the 2 month light curves. Those are very costly (in terms of CPU and disk space), and we have shown in the 4FGL paper that the 1 yr light curves capture most of the variability information. Monthly light curves, including indications of spectral variability, are now provided with similar accuracy in the Fermi-LAT Light Curve Repository.<fnref n="85" /> Detection of transient sources on monthly timescales is covered by 1FLT<fnref n="86" /> (Baldini et al. 2021). Transients on weekly timescales are covered by FAVA<fnref n="87" /> (Abdollahi et al. 2017). Phase-folded spectral analysis of periodic sources (pulsars and binaries) is provided in dedicated papers.

The condition for switching the spectral models to a curved representation is now less stringent, thereby increasing the fraction of curved spectra to about half. We improved the spectral parameterization of pulsars, quantified the residuals using a new method, and now include in the Catalog FITS file the peak energy (in $\nu F_\nu$) for each source.

The community has worked hard on the unassociated LAT sources in recent years. Optical observations allowed for classifying many blazars, and radio, optical, and γ-ray timing allowed for detecting new pulsars and binary systems. We have included all this information in the source associations. The DR2 and DR3 results provide about 900 more unassociated sources to look into. Unassociated sources close to the Galactic plane are spatially clustered and have distinct spectral characteristics. We have singled them out by a specific flag.

Section 2 describes the data and the updates to the diffuse model, Section 3 describes the updates to the analysis including the new extended sources, Section 4 describes the results including a comparison to DR1, Section 5 describes the

updates to the associations, and Section 6 discusses the unassociated sources, with particular emphasis on the category of soft clustered sources close to the Galactic plane.

## 2. Instrument and Background

### 2.1. The LAT Data

The data for the 4FGL-DR3 catalog were taken during the period 2008 August 4 (15:43 UTC) to 2020 August 2 (8:33 UTC) covering 12 years. During most of this time, Fermi was operated in sky-scanning survey mode, with the viewing direction rocking north and south of the zenith on alternate orbits such that the entire sky is observed every ~3 hr. Starting on 2018 March 16, the Fermi spacecraft was put in safe hold after one solar array drive became stuck. Scientific operations of the LAT were interrupted for more than 3 weeks, by far the longest missing time interval since 2008, and restarted on April 8 in partial sky-scanning mode.<fnref n="88" /> The Sun more rarely enters the field of view, and during some phases of the ~53-day precession period of the orbit, the entire sky is not covered every 3 hours. This hampers the solar observations and short-term light curves, but has little impact on the integrated sky coverage.

As in 4FGL, intervals around solar flares and bright gamma-ray bursts (GRBs) were excised. During the additional 4 years, 108 ks were cut due to two successive bright solar flares in 2017 September, and 2.4 ks around five new bright GRBs. The current version of the LAT data is Pass 8 P8R3 (Atwood et al. 2013; Bruel et al. 2018). For DR2 and DR3, we used the P8R3_V3 IRFs instead of P8R3_V2 used to build DR1. The relative effective areas for the point-spread function (PSF) event types were recalibrated on bright sources,<fnref n="89" /> reducing the systematics from 20% to 5%. The energy range remains 50 MeV to 1 TeV. In 4FGL, the event types were not treated separately above 10 GeV. We now decompose the events into four components up to 30 GeV, resulting in 19 components in total (Table 1), for better precision but at the cost of more disk space and CPU.

### 2.2. Model for the Diffuse Gamma-Ray Background

The interstellar emission is a very important ingredient when characterizing the point sources, because it is structured and 90% of the photons in the Galactic plane are of diffuse origin. We used exactly the same model as in 4FGL.<fnref n="90" /> The isotropic spectrum was updated to P8R3_V3. The model for the emission of the Sun and the Moon was also kept the same, neglecting the modulation of the emission along the solar cycle. It became necessary, however, to generate a different effective model for each year (for the light curves). After the solar array drive anomaly during the tenth year, the coverage of the Sun changed significantly (Figure 1). For intervals of about 2 weeks, the Sun does not enter the field of view at all.

## 3. Construction of the Catalog

Most of the steps were identical to 4FGL (Abdollahi et al. 2020). As before the extended sources are taken from the

<footnote n="85">⁸⁵ See https://fermi.gsfc.nasa.gov/ssc/data/access/lat/LightCurveRepository/.</footnote>

<footnote n="86">⁸⁶ See https://www.ssdc.asi.it/fermi_iflt.</footnote>

<footnote n="87">⁸⁷ See https://fermi.gsfc.nasa.gov/ssc/data/access/lat/FAVA/.</footnote>

<footnote n="88">⁸⁸ See https://fermi.gsfc.nasa.gov/ssc/observations/types/post_anomaly/.</footnote>

<footnote n="89">⁸⁹ See https://fermi.gsfc.nasa.gov/ssc/data/analysis/LAT_caveats_p8r3v3.html.</footnote>

<footnote n="90">⁹⁰ See https://fermi.gsfc.nasa.gov/ssc/data/access/lat/BackgroundModels.html.</footnote>
<figure><img src="image_4.png" /><figcaption>**Figure 1.** Model for the inverse Compton emission of the Sun at 1 GeV averaged over one year, in equatorial coordinates. The top and bottom rows are, respectively, years 12 (2019–2020) and 9 (2016–2017), chosen to be after and before the solar array drive anomaly. The color bar is in sqrt scale, and the intensity unit is ph cm⁻² s⁻¹ sr⁻¹ MeV⁻¹.</figcaption></figure>

<table><caption>**Table 1**
4FGL-DR3 Summed Likelihood Components</caption><thead><tr><th rowspan="2">Energy Interval (GeV)</th><th rowspan="2">NBins</th><th rowspan="2">ZMax (deg)</th><th rowspan="2">Ring Width (deg)</th><th colspan="5">Pixel Size (deg)</th></tr><tr><th>PSF0</th><th>PSF1</th><th>PSF2</th><th>PSF3</th><th>All</th></tr></thead><tbody><tr><td>0.05–0.1</td><td>6</td><td>80</td><td>7</td><td>...</td><td>...</td><td>...</td><td>0.6</td><td>...</td></tr><tr><td>0.1–0.3</td><td>5</td><td>90</td><td>7</td><td>...</td><td>...</td><td>0.6</td><td>0.6</td><td>...</td></tr><tr><td>0.3–1</td><td>6</td><td>100</td><td>5</td><td>...</td><td>0.4</td><td>0.3</td><td>0.2</td><td>...</td></tr><tr><td>1–3</td><td>5</td><td>105</td><td>4</td><td>0.4</td><td>0.15</td><td>0.1</td><td>0.1</td><td>...</td></tr><tr><td>3–10</td><td>6</td><td>105</td><td>3</td><td>0.25</td><td>0.1</td><td>0.05</td><td>0.04</td><td>...</td></tr><tr><td>10–30</td><td>5</td><td>105</td><td>2</td><td>0.15</td><td>0.06</td><td>0.04</td><td>0.03</td><td>...</td></tr><tr><td>30–1000</td><td>8</td><td>105</td><td>1</td><td>...</td><td>...</td><td>...</td><td>...</td><td>0.03</td></tr></tbody></table>

**Note.** We used 19 components (all in binned mode) in the summed-likelihood approach (Section 3.3) for DR3. Components in a given energy interval share the same number of energy bins, the same zenith angle selection, and the same region of interest (RoI) core, but have different pixel sizes in order to adapt to the PSF width. Each filled entry under Pixel Size corresponds to one component of the summed log-likelihood. NBins is the number of energy bins in the interval, ZMax is the zenith angle cut, Ring Width refers to the difference between the RoI core and the extraction region (see Abdollahi et al. 2020, for details).

literature. We only search for point sources. We use the same Test Statistic (TS), $TS = 2 \ln(\mathcal{L}/\mathcal{L}_0)$, to quantify how significantly a source emerges from the background, comparing the maximum value of the likelihood function $\mathcal{L}$ including the source in the model with $\mathcal{L}_0$, the value without the source. The 1 yr light curves were generated in exactly the same way.

### 3.1. Extended Sources

Several extended sources (existing and new) were the subject of publications since the 4FGL paper. They were not considered in DR2, but were considered for DR3:

1. We adopted the Gaussian model proposed by Devin et al. (2020) for the broad supernova remnant (SNR) G150.3 +4.5, after we confirmed that it is a better fit than the disk model used in DR1.

2. Tibaldo et al. (2018) have studied the Vela X pulsar wind nebula (PWN) and reported that the radio template was superior to the geometric templates. We found that it does not provide a better fit than the DR1 model, but it eliminates six dubious point sources, so it is simpler and we adopted it. The same authors reported a spatial evolution as a function of energy, but we found that their
<table><caption>**Table 2**
Extended Sources Entered or Modified in the 4FGL-DR3 Analysis</caption><thead><tr><th>DR3 Name</th><th>Extended Source</th><th>Origin</th><th>Spatial Form</th><th>Extent (deg)</th><th>Reference</th></tr></thead><tbody><tr><td>J0425.6+5522e</td><td>SNR G150.3+4.5</td><td>3FHL</td><td>Gaussian</td><td>1.36</td><td>Devin et al. (2020)</td></tr><tr><td>J0834.3-4542e</td><td>Vela X</td><td>2FGL</td><td>Map</td><td>1.00</td><td>Tibaldo et al. (2018)</td></tr><tr><td>J1000.0-5312e</td><td>SNR G279.0+1.1</td><td>New</td><td>Ring</td><td>0.30, 1.44</td><td>Araya (2020)</td></tr><tr><td>J1640.7-4631e</td><td>HESS J1640-465</td><td>New</td><td>Gaussian</td><td>0.11</td><td>Mares et al. (2021)</td></tr><tr><td>J1824.4-1350e</td><td>HESS J1825-137</td><td>2FGL</td><td>Gaussian</td><td>0.85</td><td>Grondin et al. (2011)</td></tr><tr><td>J2044.9+5029e</td><td>HB 21</td><td>3FGL</td><td>Disk</td><td>0.83</td><td>Ambrogi et al. (2019)</td></tr><tr><td>J2226.7+6052e</td><td>SNR G106.3+2.7</td><td>New</td><td>Disk</td><td>0.25</td><td>Xin et al. (2019)</td></tr></tbody></table>

**Note.** List of updated and new sources that have been modeled as spatially extended. The Origin column gives the name of the Fermi-LAT catalog in which that extended source was first introduced (with a different template). The Extent column indicates the radius for Disk (flat disk) sources, the 68% containment radius for Gaussian sources, the inner and outer radii for Ring (flat annulus) sources, and an approximate radius for Map (external template) sources.

high-energy disk is not significant when using the radio template.

3. Araya (2020) has shown that the cluster of DR1 sources around the SNR G279.0+1.1 is actually an extended LAT source. We adopted their ring model.

4. HESS J1640-465 is seen as extended by H.E.S.S., but was considered a point source in the DR1 catalog. We adopted the H.E.S.S. Gaussian template, after it was shown by Mares et al. (2021) that it is a better fit to the LAT data as well.

5. Two papers discussed the Fermi-LAT emission from the PWN HESS J1825-137. The morphology considered by Araya et al. (2019) was too complex for our catalog. We tried the broader Gaussian proposed by Principe et al. (2020) but found it worsened the fit. We did change the template, however, because we noticed that the one that we used since 2FGL was not the final version of Grondin et al. (2011).

6. Ambrogi et al. (2019) have shown that the previous template for the SNR HB 21 was too broad, because it mistakenly included a separate point source. We adopted their smaller template.

7. Xin et al. (2019) reported the detection with Fermi-LAT of the extended VERITAS source VER J2227+608, corresponding to the SNR G106.3+2.7 next to the much brighter pulsar PSR J2229+6114. We have included it in the DR3 model.

Table 2 lists the source name, origin, spatial template, and the reference for the dedicated analysis. These sources are tabulated with the point sources, with the only distinctions being that no position uncertainties are reported and their names end in e (see the Appendix). The new and modified extended sources have DataRelease set to 3. Unidentified point sources inside extended ones are indicated as “xxx field” in the ASSOC2 column of the catalog, where “xxx” is the name of the extended source. A total of 78 extended sources are considered in DR3.

## 3.2. Detection and Localization

The source detection procedure followed the same lines as in 4FGL. It used pointlike (Kerr 2010) and a specialized diffuse model in which the nontemplate features are estimated differently. It started from the previous list of sources (DR1 for DR2, DR2 for DR3), relocalized them over 12 yr of data, looked for peaks in the residual TS maps generated for several spectral shapes, introduced those in the model, refit, and

iterated over the full procedure. The resulting list (called uw1216) contained 11,128 seed sources at TS > 10.

That procedure naturally resulted in changing the positions of all sources. Since we wanted the catalog to be incremental, we forced most DR1 and DR2 sources back to their original positions (consistent with their names). This requires associating the uw1216 seed sources with the DR2 catalog. To do that, we started by defining pairs in which each member is the nearest counterpart of the other, within at most 1°. We repeated that procedure after removing the pairs found at the first step, in order to handle faint sources close to much brighter ones (in which the nearest counterpart can be the bright source). We rejected 58 pairs in which the distance between the two members was larger than the largest of the two 99.9% position errors (approximately 1.5 times the 95% error) and 0°.15. We eliminated the uw1216 seeds in a pair with a DR2 source, and formed the input list from the DR2 sources (including 139 with no counterpart in uw1216) and the remaining uw1216 seeds.

In DR2 we actually used the original positions of all DR1 sources strictly, and kept all of them. However a number of sources became much brighter in recent years than they were during the first 8 yr, so they are now much better localized. In recognition of this, we adopted for DR3 the new error ellipse when its area was at least 10 times smaller and TS at least 10 times larger than the same quantities in DR1. We used those very strict conditions to preserve the incremental nature of the catalog. This selected 10 sources, clear outliers when looking at the distribution of ratios of ellipse areas and TS (Figure 2). All are associated with blazars (very variable sources). One of them was in uw8606 (list of seeds for DR1) but too faint to be included in DR1 (it appeared in DR2), so only nine DR1 sources were updated. They have a 4FGL counterpart in the ASSOC_FGL column, but have DataRelease set to 3 because they have moved.

We have also adopted the new error ellipses of one DR1 source close to the modified extended source HESS J1825-137 and three faint soft sources (two from DR1, one from DR2) whose 95% semimajor axes $R_{95}$ were larger than 1° in DR1 or DR2. Three DR1 sources with broad error ellipses $R_{95} \sim 0°.4$ were split into two better localized seeds inside the original error ellipse. We replaced the DR1 sources by the two seeds. The same applied to two DR2 sources. Conversely we replaced two DR1 sources (4FGL J1831.5-0935 and J1830.2-1005) by a single uw1216 seed compatible with PSR J1831-0952.

We applied the same procedure as in 4FGL to eliminate seeds too close to a bright source and inside extended sources. Since we changed the extended sources (Section 3.1), this
<figure><img src="image_5.png" /><figcaption>**Figure 2.** Comparison of the *pointlike* localizations of uw1216 seeds (for DR3) with those of uw8606 (for DR1). The *abscissa* is the ratio of TS values, whereas the *ordinate* is the inverse of the ratio of the areas of the error ellipses (such that larger is better). Sources in the upper right corner (materialized by the dashed lines) are stronger and much better localized in uw1216, and their positions and error ellipses were updated in DR3. The source at TS ratio around 300 is the bright γ-ray binary LS 5039, which appeared double in uw8606 because of its complex spectrum (resulting in low TS) but was well localized.</figcaption></figure>

resulted in deleting several DR1 sources as well, not only new seeds. More precisely, 4FGL J1640.6–4632 was replaced by HESS J1640–465, five DR1 and two DR2 sources were deleted in or close to Vela X, three DR1 and one DR2 sources in the new SNR G279.0+1.1. We also deleted the faint 4FGL J1814.1–1710 next to the extended source HESS J1813–178 (even though this one did not change) because its TS decreased from 30 in DR1 to 16 in DR3, and it was not spectrally distinct from HESS J1813–178.

Following the *catXcheck* verification (Section 3.6), we have applied several modifications to fit the data better. We replaced 4FGL J0857.7–4507 (in one of the clusters described in Section 6.2) by a stronger uw1216 source about 0°.5 away. We replaced 4FGL J1750.0–3849 just outside the extended source FHES J1741.6–3917 by a stronger point source just inside it.

We have looked at individual faint (TS < 25) sources with a very soft (Γ > 3.1) spectrum, a large error ellipse (R$_{95}$ > 0°.5), or a stronger close neighbor in uw1216 and deleted three DR1 (4FGL J0533.9+2838, J0313.6–7508, and J2326.5+8555) and three DR2 (4FGL J0750.0+7140, J1112.0+1021, and J1752.2–3002) sources. Nearby seeds appeared in DR3 and fit the data better. Similarly, we have deleted a faint very hard (Γ < 1.3) DR2 source (4FGL J1830.3–1601) that was superseded by a brighter nearby DR3 seed.

As before we manually added 11 known LAT pulsars that could not be localized by the automatic procedure without phase selection. Only one of those (PSR J1731–4744) reached TS > 25 in DR3. It appears with a NULL error ellipse. The 10 others were discarded. A total of 5029 new seeds were entered to the *gtlike* source characterization for DR3 in addition to the 5690 remaining DR1 and DR2 point sources.

We reassessed the systematic corrections to localization. Up to 4FGL, they were defined globally over the entire sky, so their values were dominated by high-latitude sources (92% of associated point sources in DR1 are at |b| > 5°). However, confusion and strong interstellar emission make localization much more difficult in the Galactic plane. In DR2 we estimated localization systematics separately at |b| < 5°, based on 110 pulsars and low-latitude active galactic nuclei (AGN). We found that they should indeed be increased there, to 27" for the

<table><caption>**Table 3**
4FGL-DR3 Priors on Diffuse Parameters</caption><thead><tr><td>Parameter Name</td><td>Mean</td><td>Std Dev</td></tr></thead><tbody><tr><td>Galactic norm at 1 GeV, K</td><td>0.97</td><td>0.03</td></tr><tr><td>Galactic spectral bias, Γ<sub>d</sub></td><td>-0.01</td><td>0.02</td></tr><tr><td>Isotropic normalization</td><td>1.00</td><td>0.10</td></tr></tbody></table>

**Note.** The Galactic interstellar emission model is modulated by a power law K(E/E$_{0}$)$^{-\Gamma_{d}}$, in which E$_{0}$ is set to 1 GeV. K tends to be smaller than 1. This is probably because of the sources that we added since the model was devised (just before 4FGL).

absolute 95% error, and 1.37 for the systematic factor. We did not change the systematics in the high-latitude sky (25" and 1.06). The DR2 systematics are used again for DR3.

### 3.3. Thresholding and Light Curves

We used the Fermi Tools 1.4.7 analysis suite.$^{91}$ Compared to the Science Tools v11r7p0 that were used in 4FGL, the main difference is the new *edisp bins* setting for energy dispersion, which allows for broadening of the energy interval considered in the model space. The Science Tools setting was equivalent to *edisp bins* = -1, not enough to fully fill the first and last data bins for the small energy bins that we use (10 per decade). We used *edisp bins* = -2 instead, considering two bins for energy dispersion on each side beyond the data bins. The acceptance for the PSF3 event type increases by a factor 2.4 between 50 and 100 MeV, so to account for this rapid change, we used six energy bins in the first component (Table 1) instead of three and set *edisp bins* to -4 accordingly.

The likelihood weights were recomputed over 12 yr of data, resulting in slightly smaller weights everywhere. We reoptimized all RoIs, resulting in 1988 RoIs containing up to 10 sources in their core. The maximum of eight sources in the RoI core enforced in 4FGL was not optimal (it led to too many small neighboring RoIs and too many iterations).

For DR1 and DR2, the normalization of one of the diffuse components sometimes became very large (close to the maximum, set to 2) or very small (close to 0) in regions where that component contributed little. This was particularly obvious for the isotropic component in the Galactic ridge. This could happen because of missing features in the interstellar emission model (IEM), leading to a higher isotropic component in the fit. This was not satisfying since the isotropic emission is supposed to be approximately constant over the sky, and the stuck parameters caused trouble in the error estimates. In order to stabilize the diffuse parameters, we replaced the hard limits with Bayesian priors. We chose the priors from the distributions of the parameters (in DR2) in regions where they are dominant (|b| < 10° for the Galactic parameters and |b| > 30° for the isotropic normalization). This led to the values reported in Table 3. The intrinsic reasons for the large or small values did not go away, of course, so the resulting isotropic normalizations were still larger than one in the Galactic ridge, but they rarely exceeded 1.5. The worst RoIs around (l, b) ~ (315°, 3°) reached 1.7.

A major difference with the 4FGL procedure is that we retained all previous 4FGL sources (DR1 and DR2, except

<footnote>$^{91}$ See https://fermi.gsfc.nasa.gov/ssc/data/analysis/documentation/.</footnote>
those discussed in Section 3.2), so that they would not be deleted from the model even if they have TS < 25. The resulting catalog contains 6659 entries, among which 1610 are new (1607 point sources plus the three extended sources discussed in Section 3.1), and 181 have TS < 25 (112 DR1 and 69 DR2 sources).

The light curves over 1 yr bins were recomputed (even the first 8 yr because adding new sources changes the fluxes of DR1 sources). The procedure was very similar to 4FGL, except we considered events up to 1 TeV, and we used a different model for the Sun and Moon at each year (see Section 2.2). Neither change has a major effect. Recomputing light curves over 2 month bins would have been very time consuming, and we showed in the 4FGL paper that the finer time binning was not critical to detecting variability. Monthly (and shorter intervals) light curves are provided in the Fermi-LAT Light Curve Repository.

The threshold on Variability_Index is now 24.725 (same 1% probability of false detection, applied to 12 bins). The number of significantly variable sources increased from 1443 to 1695 (among which 158 are DR2 or DR3 sources). The fraction of variable sources decreased because part of the variability information in 4FGL came from the 2 month light curves. The distribution of fractional variability (Frac_Variability in the FITS file) is similar to 4FGL, peaking between 50% and 90%. The blazars all have fractional variability larger than 10% except 4C +55.17 (4FGL J0957.6 +5523). The flux of that source, first discussed by McConville et al. (2011), decreased by about 13% over 10 yr, not showing any flare.

### 3.4. Spectral Fitting and Energy Distributions

The spectral analysis proceeded along the same lines as in 4FGL, with five important differences detailed below: we changed the parameterization for pulsars, we fit more sources with four parameters, we changed the TS\_curv threshold for using a curved spectral shape, we added new columns reporting the peak energy in $\nu F_{\nu}$, and we added one spectral bin to the spectral energy distributions (SEDs).

Most significantly curved sources are still fit with a lognormal function (LogParabola under SpectrumType in the FITS table, hereafter LP):

$$ \frac{dN}{dE} = K \left( \frac{E}{E_0} \right)^{-\alpha - \beta \ln(E/E_0)} \quad (1) $$

The significantly curved pulsars are still fit with a subexponentially cutoff power law (PLSuperExpCutoff under SpectrumType in the FITS table, hereafter PLEC). However we noticed that the parameters of the parameterization used in 4FGL (PLSuperExpCutoff2 in the Fermi Tools) were strongly correlated, particularly when the exponential index ($b$ in Equation (2)) was free. This hampered convergence and precluded freeing $b$ except in the very brightest sources. To work around this, we devised a new parameterization (PLSuperExpCutoff4 in the Fermi Tools):

$$ \frac{dN}{dE} = K \left( \frac{E}{E_0} \right)^{\frac{d}{b} - \Gamma_S} \exp \left[ \frac{d}{b^2} \left( 1 - \left( \frac{E}{E_0} \right)^b \right) \right] \quad (2) $$

<figure><img src="image_6.png" /><figcaption>Figure 3. Ratio of the error on the $\Gamma_S$ PLEC parameter (Equation (2)) in DR3 to that on the asymptotic low-energy index $\Gamma = \Gamma_S - d/b$ in DR2 (diamonds). The dots show the same error ratio on the $\alpha$ LP parameter (Equation (1)), to gauge the effect of increased exposure (small). The large scatter at significance larger than 100 is due to the fact that those sources have free exponential index $b$ in DR3.</figcaption></figure>

$$ \frac{dN}{dE} = K \left( \frac{E}{E_0} \right)^{-\Gamma_S - \frac{d}{2} \ln \frac{E}{E_0} - \frac{db}{6} \ln^2 \frac{E}{E_0} - \frac{db^2}{24} \ln^3 \frac{E}{E_0}} \quad \text{for} \\ \left| b \ln \frac{E}{E_0} \right| < 10^{-2}, \quad (3) $$

in which the normalization $K$ is directly the flux density at the reference energy $E_0$, and the shape parameters are the spectral slope $\Gamma_S = d \log(dN/dE)/d \log E$ and the spectral curvature $d = d^2 \log(dN/dE)/d(\log E)^2$ at $E_0$.

The development in Equation (3) shows explicitly that for $b=0$ the expression converges to LP with $d=2\beta$. For $b<0$ the shape (in log) is reversed (it decreases exponentially toward low energies and as a PL toward high energies). The parameters $K$, $\Gamma_S$, $d$, and $b$ appear as PLEC_Flux_Density, PLEC_IndexS, PLEC_ExpfactorS, and PLEC_ExpIndex in the FITS table, respectively. The reference energy $E_0$ is chosen so that the error on $K$ is minimal and appears as Pivot_Energy in the FITS table. This parameterization proved indeed much more stable, and the correlation between parameters is considerably reduced.

Figure 3 shows that on average $\Gamma_S$ is much better defined (error ratio below 1) than the low-energy asymptotic index $\Gamma = \Gamma_S - d/b$ used in 4FGL. On average the error on $\alpha$ ($\Delta\alpha$) decreased by only 2.3% between 10 and 12 yr, whereas $\Delta\Gamma_S$ is on average 0.626 $\Delta\Gamma$. Note that 15 pulsars had no $\Delta\Gamma$ at all in DR2 because $\Gamma$ was stuck at 0. This problem does not exist in the new parameterization. The small error ratio at significance below 40 corresponds to PSR J0622+3749. It had $\Gamma=0.05$ in DR2, very close to the boundary. At low significance, a fraction of the pulsars have $\Delta\Gamma$ in DR2 that is similar or less than $\Delta\Gamma_S$ in DR3. The DR2 errors may have been underestimated because of the large correlation between parameters.

Since the new PLEC parameters describe the same mathematical quantities (spectral slope and curvature) as the LP parameters (except the curvature is $2\beta$), it is interesting to compare them. Both spectral models are fit to all sources, so this information is readily available. Figure 4 shows that, as expected, the correlation between $\alpha$ and $\Gamma_S$ is quite tight. The correlation between $2\beta$ and $d$ is not as tight, because the curvature has larger error. On average $d<2\beta$, because the
<figure><img src="image_7.png" /><figcaption>**Figure 6.** Ratio between the energy flux (0.1 to 100 GeV) in the PL approximation and that in the LP or PLEC representation, for those sources that are considered significantly curved, outside the Galactic plane. Plus signs refer to the LP model, diamonds to the PLEC model. The scatter increases toward low significance because of the larger statistical fluctuations.</figcaption></figure>

<figure><img src="image_8.png" /><figcaption>**Figure 7.** SED of the BL Lac object PG 1218+304 (4FGL J1221.3+3010) for DR2 (purple crosses and LP fit) and DR3 (green diamonds and LP fit).</figcaption></figure>

the ridge (green) diverges toward low energies at a threshold of 9 (what was actually used in DR2). The global spectrum is better behaved when the TS<sub>curv</sub> threshold is lower.

The implication is that we should choose a TS<sub>curv</sub> threshold as low as possible while preserving global convergence of the model. For DR3 we set the threshold to 4 as a compromise because the test with a threshold at 1 required manual intervention, and global convergence could not be guaranteed. This is still a significant improvement as demonstrated in Figure 5 (right).

Thanks to the improved statistics and the lowered TS<sub>curv</sub> threshold, many more sources are considered significantly curved in DR3 than in DR1. The number of LP spectral models is now 3130 (among which 650 are new sources), while the number of PLEC models is now 258 (among which 2 are new sources). The fraction of sources represented by a curved spectral model increased from 30% to 51% overall.

Changing the spectral models for many sources from PL to curved also has an impact on their energy flux estimate. The PL flux is nearly always larger, due to the integration down to 100 MeV for soft sources (or up to 100 GeV for hard sources Γ < 2). This is quantified in Figure 6. Well-measured pulsars (diamonds at high significance) have on average a 50% larger flux estimate in the PL model than in the PLEC model. The well-measured blazars (plus signs) are in the lower branch at a ratio

around 1.1, because their spectra are much less curved than those of pulsars. The diamonds in the lower branch are the very bright blazars fit with PLEC (see above). A ratio less than one occurs when the peak of the curved spectrum is close to either 100 MeV or 100 GeV (the boundaries of the energy interval over which the energy flux is extracted), because the full analysis interval (50 MeV to 1 TeV) extends beyond those boundaries and maintains the PL fit below the peak of the curved model.

To provide more information on curved spectra, we now report systematically in the catalog the peak energy in νF<sub>ν</sub> and its uncertainty for all sources (including those not significantly curved) and both models as (Unc_)LP_EPeak and (Unc_)PLEC_EPeak:

$$E_{\text{peak}}(\text{LP}) = E_0 \exp\left(\frac{2 - \alpha}{2\beta}\right) \qquad (4)$$

$$E_{\text{peak}}(\text{PLEC}) = E_0 \left( 1 + \frac{b}{d} (2 - \Gamma_S) \right)^{1/b} \qquad (5)$$

The uncertainties are obtained using the covariance matrix. We set those columns to NULL when the peak energy is undefined (when Γ ≥ 2 in PLEC or curvature is upwards) or found too far from the LAT energy range (|ln(E<sub>peak</sub>/E<sub>0</sub>)| > 10). E<sub>peak</sub>(PLEC) is undefined in only 4 PLEC spectra, but in 516 LP and 1682 PL spectra (mostly because Γ ≥ 2). E<sub>peak</sub>(LP) is always defined in PLEC and LP spectra, but remains undefined in 850 PL spectra in which β is found negative or close to 0. The uncertainty on E<sub>peak</sub>(LP) is also set to NULL when β is stuck to 1 and has no error (166 occurrences in LP spectra, 6 in PLEC spectra, and 3 in PL spectra). Even when E<sub>peak</sub> is formally defined, it should be treated with caution when Γ (PLEC) gets close to 2 or the curvature is close to 0.

In the DR1 and DR2 catalogs, we provided SEDs over N<sub>bands</sub> = 7 energy bands. Thanks to the increased statistics after 12 yr (particularly valuable at high energy where the precision is limited by counts rather than background), we can provide them over eight energy bands in DR3, extending to 1 TeV as in the 3FHL catalog (Ajello et al. 2017). We replaced the 30–300 GeV bin by a 30–100 GeV bin and a 100 GeV–1 TeV bin. We computed those two fluxes using unbinned likelihood, as in 4FGL. In the highest-energy band, 172 sources have TS > 25. An example of such a hard source is shown in Figure 7. That source has TS<sub>curv</sub>(LP) = 8.70, below the DR2 threshold of 9 but above the DR3 threshold of 4, so it is fit with the LP model in DR3. The significant differences between 1 and 10 GeV are due to the natural variability of that strong source.

The systematic uncertainties<sup>92</sup> have improved when moving from P8R3_V2 to P8R3_V3 and are now 0.1 in bands 1 and 8 and 0.05 in bands 2–7. We compute the spectral fit quality (SFQ) as before (Equation (5) in Abdollahi et al. 2020). Figure 8 shows how that quantity is distributed. Up to DR2, we flagged (Flag 10 in Table 4) sources with SFQ above the 99.9% quantile (10<sup>-3</sup> p-value) of the χ<sup>2</sup> distribution with N<sub>bands</sub> − 2 degrees of freedom (dof), considering that most sources were fit with PL spectra that have two free parameters. For DR3 that would result in 6 dof. Clearly the data do not follow the χ<sup>2</sup>(6) distribution, but are much closer to the χ<sup>2</sup>(4) distribution. This is due to two effects: first, half the sources are now fit with curved spectra (3 free parameters); second, the first (<100 MeV)

<footnote n="92">92 See https://fermi.gsfc.nasa.gov/ssc/data/analysis/LAT_caveats.html.</footnote>
<figure><img src="image_9.png" /><figcaption>**Figure 8.** Distribution of spectral fit quality over all DR3 sources (histogram). The $\chi^2$ distributions for 4 (dashed line, good fit) and 6 (solid line, bad fit) degrees of freedom are overlaid. The vertical lines correspond to a $p$-value of $10^{-3}$. The threshold for Flag 10 in Table 4 corresponds to the dashed vertical line (for 4 dof).</figcaption></figure>

and last (>100 GeV) bands do not behave as a Gaussian distribution (very large errors in the first band, and very low counts in the last band), so they contribute little to SFQ. For DR3 we flagged sources with SFQ > 18.47, the $10^{-3} p$-value of the $\chi^2(4)$ distribution. The number of sources that trigger this remains small (only 16). The largest SFQ is obtained in PSR J0205+6449 at 49.6. That pulsar indeed has a known PWN (3C 58), which is very significant above the PLEC fit in the two bands above 30 GeV. The other sources flagged with Flag 10 have one spectral band off by more than $3\sigma$.

### 3.5. Analysis Flags

As before, we attempt to flag the sources affected by systematic errors larger than the statistical errors that we report, and which should therefore be treated with more caution. This can impact all source parameters, and even call into question the existence of the faintest flagged sources with significance less than $6\sigma$ or $7\sigma$ (TS < 50). Because the largest source of systematics is the underlying diffuse emission, the sources close to the Galactic plane, where diffuse emission is both the strongest and most structured, are most affected.

The flags (Flags in the FITS file) are recalled in Table 4 (see the 4FGL paper for the detailed definitions), together with the numbers of sources flagged for each reason and their evolution since DR1. The effect of the underlying IEM was estimated by launching the procedure described in Section 3.3 for a second time using the same seeds but the previous IEM (gll_iem_v06).

We note six changes with respect to the 4FGL procedure:

1. For DR1, Flag 1 was applied to the comparison with both *pointlike* and the alternative IEM. For DR3, the comparison to *pointlike* was performed with no TS selection in order to apply Flag 3 to sources at TS < 25. Flag 1 was determined only by the comparison with the old IEM, after selecting TS > 25 on both sides.

2. Flag 2 mostly relies on the comparison between the new uw1216 *pointlike* positions of the 4FGL sources (Section 3.2) and the 4FGL positions (uw8606), because no localization was performed over 12 yr with the old IEM. The comparison with the old IEM also flags a few sources (although the seed positions are exactly the same) when among two nearby seeds one survives in DR3 and the other one with the old IEM.

<table><caption>**Table 4**
Comparison of the Numbers of Flagged Sources between DR1 and DR3, Separately for the DR1 Sources (“DR1 in DR3” Column) and the Sources Introduced in DR2 and DR3 (“Post-DR1” Column)</caption><thead><tr><th>Flag<sup>a</sup></th><th>DR1</th><th>DR1 in DR3</th><th>Post-DR1</th><th>Meaning</th></tr></thead><tbody><tr><td>1</td><td>215</td><td>151</td><td>117</td><td>TS &lt; 25 with other model or analysis</td></tr><tr><td>2</td><td>215</td><td>393</td><td>50</td><td>Moved beyond 95% error ellipse</td></tr><tr><td>3</td><td>342</td><td>352</td><td>139</td><td>Flux changed with other model or analysis</td></tr><tr><td>4</td><td>212</td><td>233</td><td>227</td><td>Source/background ratio &lt;10%</td></tr><tr><td>5</td><td>398</td><td>418</td><td>259</td><td>Confused</td></tr><tr><td>6</td><td>92</td><td>161</td><td>156</td><td>Interstellar gas clump (c sources)</td></tr><tr><td>9</td><td>136</td><td>97</td><td>71</td><td>Localization flag from <i>pointlike</i></td></tr><tr><td>10</td><td>27</td><td>44</td><td>2</td><td>Bad spectral fit quality</td></tr><tr><td>12</td><td>103</td><td>94</td><td>90</td><td>Highly curved spectrum</td></tr><tr><td>13</td><td>...</td><td>112</td><td>69</td><td>TS &lt; 25 at 12 yr</td></tr><tr><td>14</td><td>...</td><td>368</td><td>181</td><td>Soft Galactic Unassociated (Section 6.2)</td></tr><tr><td>All</td><td>1163</td><td>1433</td><td>720</td><td>Any flag (Flags &gt;0)</td></tr></tbody></table>

**Note.**

<sup>a</sup> In the FITS file, the values are encoded as individual bits in the *Flags* column, with Flag *n* having value $2^{(n-1)}$.

3. All sources have been visually screened for potential deficiencies in the underlying diffuse model. The visual screening for diffuse features (Flag 6) depends on the source parameters. It was entirely redone for DR2. To be conservative, we left that flag (and the *c* suffix) in all the sources that had it in DR1, even if they would not have been flagged in DR2, and added the flag to all the new cases, including the sources already in DR1. For DR3 we checked only the new sources and did not change Flag 6 in the DR1 and DR2 sources. Two sources in the Large Magellanic Cloud (4FGL J0517.9-6930c and J0535.7-6604c) have that flag set as well.

4. Flag 12 still applies to PLEC sources with $\Gamma = \Gamma_s - d/b \le 0$ (as well as to LP sources with $\beta$ fixed to 1, as before), even though the new PLEC parameterization (Equation (2)) has no explicit limit at $\Gamma = 0$.

5. We introduced for DR2 a new Flag (13) to flag explicitly sources from a previous version with TS < 25 in the current version. It now applies to DR2 sources as well (69 have TS < 25 in DR3).

6. For DR3, we introduce a new Flag (14) to flag explicitly soft Galactic unassociated (SGU) sources in regions with high source density (Section 6.2). The three strongest SGUs (4FGL J1702.7-5655 at $b = -9.2$, 4FGL J1823.3-1340 inside HESS J1825-137, and 4FGL J1839.4-0553 inside HESS J1841-055), all at TS > 800, are very clear point sources (good pulsar candidates) and are not flagged.

Overall the fraction of sources with any flag set increased from 23% in DR1 to 32% in DR3. Among the new DR2 and DR3 sources, 75% are flagged in the Galactic plane ($|b| < 10^\circ$), but only 26% are flagged above the Galactic plane.

### 3.6. All-sky Verification

In order to verify that the 4FGL-DR3 catalog provides a good representation of the $\gamma$-ray sky, we have performed an independent all-sky analysis, named *catXcheck*, using the
p-value statistic (PS) data/model deviation estimator developed by Bruel (2021). *catXcheck* consists of the analysis of 540 RoIs (12 × 12°, with a pixel size of 0°.1) covering the whole sky. The RoI centers lie on 19 Galactic parallels, whose latitudes range from −90° to 90° with a 9° step. The longitude step is 9° for |*b*| ≤ 18°, 10° for |*b*| = 27° and 36°, and it is 12°, 15°, 18°, 24°, and 45° for |*b*| = 45°, 54°, 63°, 72°, and 81°, respectively. The difference with respect to the catalog main analysis is that we combine all PSF event types together, using all events whose zenith angle is less than 90° to suppress contamination from Earth limb γ rays. The energy range is [100 MeV, 1 TeV], and we use 10 bins per decade.

The source model of each RoI comprises the Galactic diffuse emission, the isotropic template, the Sun and Moon steady emission templates, as well as all point-like and extended sources from DR3, within $5 + 0.015(\sigma_{\text{src}} - 4)$ degrees of the RoI border, where $\sigma_{\text{src}}$ is the source significance (Signif_Avg in the FITS version of the catalog). That formula selects sources outside of the RoI up to 5° for the faintest ones and to 15° for the brightest ones, corresponding respectively to the 68% and 95% containment radii at 100 MeV.

Since we want to check the goodness of the DR3 information, the spectral parameters of the DR3 sources are fixed to their catalog values, and the only free parameters of the binned-likelihood fit are the ones of the Galactic diffuse emission (power-law correction) and isotropic (normalization) templates. After the fit, we compute the 3D map of the number of photons predicted by the source model and compute the PS map in order to detect any significant deviation between the observed and predicted 3D maps.

As explained in Bruel (2021), the PS estimator is computed for each spatial pixel of a RoI, and rather than comparing the observed and predicted count spectra at that pixel, we compare spatially integrated count spectra around that pixel, with an integration radius varying with energy as $p_0(E/100 \text{ MeV})^{-p_1} \oplus p_2$, where ⊕ is the sum in quadrature. We use the point-source optimized integration parameterization (which is very close to the PSF 68% containment energy dependence): $p_0 = 4^\circ$, $p_1 = 0.9$ and $p_2 = 0^\circ.1$, and we merge the energy bins to have an energy binning of 0.3 in $\log_{10} E$. PS is defined such that $|\text{PS}| = -\log_{10}(p - \text{value})$, where the p-value is the probability that the statistical fluctuations can reach a level of deviation as large as the one observed in the data, under the assumption that the model represents the data. The PS sign is the sign of the sum of the count spectra residuals in units of $\sigma$. The likelihood weights are taken into account in the PS computation. We chose to use the PS map to investigate data/model deviations because, contrary to the TS map, it is sensitive to both positive and negative deviations.

For one pixel, the $3\sigma$, $4\sigma$, and $5\sigma$ thresholds correspond to $PS = 2.57, 4.20$ and $6.24$. Since the level of correlation between pixels is negligible with $p_2 = 0^\circ.1$, we only have to take into account the number of trials when deriving the all-sky significance thresholds. The number of RoIs and their centers have been chosen so that the overlap between neighboring RoIs ensures that, for any given direction in the sky, there is at least one RoI containing this direction and providing a distance to the RoI border larger than ~1°.5. As a consequence, when combining the results of the 540 PS maps, we loop over the RoIs and the RoI pixels, and we only consider the pixels for which the RoI provides the largest distance to the RoI border. This allows us to avoid double counting the sky directions. The total number of considered pixels is 4,175,468 and adding its

<figure><img src="image_10.png" /><figcaption>**Figure 9.** PS distribution for the three Galactic-latitude selections. The thin solid histogram corresponds to the expected $10^{-|x|}$ distribution. All histograms are normalized such that their integral is 1.</figcaption></figure>

<figure><img src="image_11.png" /><figcaption>**Figure 10.** All-sky PS map with “source” settings in Galactic coordinates (Mollweide projection). The color scale corresponds to $-4\sigma$ to $4\sigma$ deviations (accounting for the number of pixels in the map).</figcaption></figure>

logarithm to the one-pixel significance thresholds gives the all-sky $3\sigma$, $4\sigma$, and $5\sigma$ thresholds: $PS = 9.19, 10.82$, and $12.86$.

The PS distributions for samples of sources in three different ranges of Galactic latitude are shown in Figure 9. The high-latitude ($|b| > 35^\circ$) distribution follows the $10^{-|x|}$ expectation, whereas some deviations are visible in the low-latitude ($|b| < 5^\circ$) and mid-latitude ($5^\circ < |b| < 35^\circ$) distributions. To show the locations of these deviations, we build an all-sky HEALPix (Górski et al. 2005) map in Galactic coordinates with $N_{\text{side}} = 256$ from the 540 PS maps. Figure 10 shows the resulting all-sky PS map. Overall the data/model agreement is good over the whole sky. Using a $4\sigma$ threshold, we find two positive and five negative deviations. The positive ones are both located inside bright extended sources: 4FGL J1923.2+1408e (W 51C) and 4FGL J0822.1–4253e (Puppis A). They occur because the geometric templates that we use for those sources are far from perfect, and we do not allow additional point sources inside them (Section 3.2). Regarding the negative deviations, the one at a distance of ~1° from the Galactic center corresponds to only one pixel. The other four are clusters of several pixels with $PS < -10.82$ and correspond to extended deviations explaining the negative tail in the mid-latitude histogram in Figure 9. They all are close to large molecular clouds: one can be associated with Cepheus, one with Perseus, and two with Orion B.

Since these extended negative deviations seem to be related to the Galactic diffuse emission, we also compute the 540 PS maps with the set of PS spatial integration parameters
<figure><img src="image_12.png" /><figcaption>**Figure 11.** All-sky PS map with “diffuse” settings in Galactic coordinates (Mollweide projection). The color scale corresponds to $-4\sigma$ to $4\sigma$ deviations (accounting for the number of independent pixels in the map).</figcaption></figure>

optimized for large extended deviations ($p_0 = 5^\circ$, $p_1 = 0^\circ.8$, and $p_2 = 1^\circ$). With $p_2 = 1^\circ$, the level of correlation between pixels is expected to increase. Assuming that pixels within $0^\circ.5$ of each other are highly correlated, we divide the number of trials by 69, and the $3\sigma$, $4\sigma$, $5\sigma$ thresholds become PS = 7.4, 9, 11. The corresponding all-sky PS map, shown in Figure 11, exhibits many deviations above $4\sigma$. A preliminary analysis of the negative ones shows that many of them are correlated with some components of the IEM, especially the local CO template (Galactocentric radius range 7–9 kpc), which may indicate that allowance for variation of the $X_{CO}$ factor among the local clouds is necessary. The positive deviations are less numerous and less significant. Like the negative ones, some are likely due to imperfections in the modeling of the Galactic diffuse emission, but some could correspond to new extended sources.

Although a complete analysis of these extended deviations is outside the scope of this paper, we try to quantify their impact on the DR3 catalog. In order to do so, we look at the 540 PS maps and count the number of pixels with $|PS| > 7.4$ and derive the expected number of DR3 sources in these regions from the average source density in the $|PS| < 7.4$ region. In the 40 RoIs centered along the Galactic plane, the fraction of negative and positive deviation pixels is 2% and 0.9%, respectively. The numbers of observed and expected DR3 sources are 15 and 24.2 in the negative deviation pixels, and 33 and 10.1 in the positive deviation pixels. In the RoIs farther from the Galactic plane, the fraction of negative and positive deviation pixels is 0.15% and 0.03%, respectively. There are no DR3 sources in the negative deviation pixels, while 8.4 are expected. We conclude that the catalog results are biased in the negative and positive deviation regions, but the effect is small and the number of missed or spurious sources is at most $\sim 20$.

## 4. The 4FGL-DR3 Catalog

The catalog is available online,<fnref n="93" /> together with associated products. It contains 6659 entries (1607 new since DR1). The source designation is 4FGL JHHMM.m+DDMM. The format is similar to 4FGL (except it does not contain the 2 month light curves) but not identical. It is detailed in the Appendix. A new column (*DataRelease*) is set to 1 for the unchanged DR1 sources, 2 for DR2 sources, and 3 for DR3 sources. Sources existing in DR1 or DR2 but whose position (and therefore name) changed (Section 3.2) have *DataRelease* set to 3.

<figure><img src="image_13.png" /><figcaption>**Figure 12.** Distributions of the 4FGL-DR3 Test Statistic at high latitude ($|b| > 10^\circ$) for the sources already in DR1 (dashed) and those new to DR2 or DR3 (solid). The vertical dotted line is the threshold at TS = 25. The solid histogram extends below the threshold because 69 DR2 sources have TS < 25 in DR3, including 46 at $|b| > 10^\circ$.</figcaption></figure>

The DR1 or DR2 names of those sources are given in the ASSOC_4FGL column of the FITS file. They are considered DR1 or DR2 sources in the text, tables, and plots.

### 4.1. The 4FGL-DR1 Sources

The TS of the 5048 remaining DR1 sources in DR3 increased on average compared to DR1. Among sources that are fit with the PL model in both catalogs (no TS boost from switching to a curved model), TS increased on average by 33% at high latitude ($|b| > 10^\circ$) and by 22% closer to the Galactic plane. The difference between the two is due to the lower log-likelihood weights, which partly offset the exposure increase close to the Galactic plane where diffuse emission is strong, particularly for soft sources. Even at high latitude, the TS ratio is less than expected from the exposure increase (50%). This is not due to the lower weights (the TS ratio does not depend on index) but to the selection bias (sources tend to be brighter in the interval in which they were defined) and signal-splitting with new sources. The statistical uncertainties on the parameters decreased by 15% at high latitude and 9% close to the plane.

The random character of adding new data and the variable nature of the γ-ray sources inevitably lead to broadening the TS distribution, in addition to shifting it to larger values. The sharp selection threshold at TS = 25 in DR1 becomes a tail of sources at TS < 25 in DR3. Figure 12 shows the TS distribution outside the Galactic plane, where 74 DR1 sources have TS < 25. Over the entire sky, 112 DR1 sources have TS < 25 in DR3, among which 15 have TS < 16, and the smallest TS is 9.8 (4FGL J1816.4–0057).

The energy flux decreased on average by 5% ($0.2\sigma$) between DR1 and DR3 (considering only PL sources to avoid the effect of switching to a curved model). This is again a selection bias dominated by faint sources: close to the detection threshold, DR1 sources can get either brighter or fainter in DR3, but faint sources in DR3 cannot be fainter in DR1 because they would have TS < 25. Real variability is apparent for bright sources (TS > 1000). The scatter on the energy flux ratio is $5\sigma$ on variable DR1 sources, but only $1\sigma$ on nonvariable ones. There is no offset and less scatter (only $0.8\sigma$ overall) on the difference of PL indices because source fluxes vary strongly, but spectral shapes only slightly.

<footnote n="93">⁹³ See https://fermi.gsfc.nasa.gov/ssc/data/access/lat/12yroprism/.</footnote>
<figure><img src="image_14.png" /><figcaption>**Figure 13.** Left: distributions of the 4FGL-DR3 energy flux at high latitude ($|b| > 10°$) for the sources already in DR1 (dotted) and those new to DR2 or DR3, separately for those considered variable (dashed) or not (solid). The energy flux is obtained with the best spectral model, which can be either a curved model or a power law (Section 3.4). Right: same distributions for the PL index (always obtained with the power-law model). The DR1 distribution was broadened to account for the larger uncertainties of the index values for the new sources, so the comparison is fair.</figcaption></figure>

<figure><img src="image_15.png" /></figure>

## 4.2. The New DR2 and DR3 Sources

The 1607 new DR2 and DR3 sources are on average very close to the detection threshold (median TS of 34). This is obvious in Figure 12. There are a few newly active bright blazars at TS > 100, and even three at TS > 800 (4FGL J0154.5−6604, J1544.3−0649, and J2250.0−1250 associated with PKS 0153−663, NVSS J154419−064913, and PKS 2247−131 respectively), but the majority are faint sources that rise above the detection threshold thanks to the additional four years. Outside the Galactic plane, the confusion limit is still distant, because ever-deeper exposures reach to higher energies at which the PSF is better. Among the 1042 faint (TS < 100) new sources at $|b| > 10°$, only 7% have the confusion flag (Flag 5 in Table 4) set.

Figure 13 (left) shows that the energy flux distribution of the new sources peaks between 1 and $2 \times 10^{-12}$ erg cm$^{-2}$ s$^{-1}$. The majority (140) of the 156 new DR2 and DR3 sources considered variable (mostly blazars) are outside the Galactic plane ($|b| > 10°$). Their average flux is larger than that of the nonvariable sources. This is largely a selection bias (variability cannot be detected in the very faint sources).

Figure 13 (right) compares the distribution of PL indices of the new sources with that of the DR1 ones. The variable new sources are soft on average, similar to flat-spectrum radio quasars (FSRQs). The same effect is observed in the 1FLT catalog (Baldini et al. 2021), even more strongly because detections in 1FLT are based on one-month data sets. The nonvariable new sources have a peak on the soft side as well (probably similar to the variable ones, but not bright enough to be significantly variable). But they also have a shoulder on the hard side. This is because adding more exposure is more beneficial at high energies (in the count-dominated regime) than at low energies, and therefore helps with detecting hard sources (like BL Lac objects).

## 4.3. Step-by-step from DR1 to DR3

To understand the effects of the analysis changes introduced in DR3 with respect to 4FGL, we have considered their effects on the DR2 data set (10 years). To that end, we started with the same seeds as the DR2 catalog, changed each element in sequence (in the order of the list below), and compared each intermediate result with the previous one.

1. We first switched from the Science Tools to the Fermi Tools, and from the P8R3_V2 to the P8R3_V3 data and Instrument Response Functions (changing the isotropic templates accordingly). This led to essentially no change at all, as expected.

2. Setting `edisp bins` to -2 (Section 3.3) when accounting for energy dispersion had a significant effect. Prior to that, the model counts were underestimated by about 5% in the first energy band (50 to 100 MeV), so the model fit was too high by the same amount at low energy. We found 16 fewer sources with the new settings, due to a slight TS reduction in soft sources. The energy flux decreased by 2% on average. The largest effect was spectral: we found 155 more LP spectra (+10%). Hard PL sources ($\Gamma < 2.1$) were unchanged, but $\Gamma$ decreased more and more for softer sources, reaching $\Delta\Gamma = -0.1$ at $\Gamma = 3$. Note that this is specific to the catalog setting that uses many components for the summed likelihood in relatively narrow energy bands (Table 1). Analyses using few broader energy bands are not as sensitive to `edisp bins`.

3. Introducing diffuse priors had very little effect on the sources. We found 7 fewer sources and 11 fewer LP spectra. In the Galactic ridge (where the effect of the priors is strongest), the average energy flux of sources decreased by 2% (0.1$\sigma$), and the average PL index decreased (became harder) by 0.005 (0.06$\sigma$). This is completely innocuous.

4. Lowering the $TS_{curv}$ threshold from 9 to 4 (Section 3.4) had an important effect on the sources, besides the collective effect illustrated in Figure 5 (right). The number of LP spectra increased by 1000 (+58%). The total number of sources increased by 127 (+2%) because more LP sources allowed more soft neighbors to reach the detection threshold. The energy fluxes of sources that switched from PL to LP decreased by 1.8$\sigma$ (well above the statistical precision) as illustrated in Figure 6. Their neighbors picked up most of that flux, so their energy fluxes increased. They also became somewhat softer in the process (by 0.02 on average). On the other hand, the neighbors of those neighbors became a little harder, so that the global average PL index did not change.
In conclusion, the change with the greatest consequences was lowering the TS<sub>curv</sub> threshold. We are confident that the new setting is closer to the physical reality.

**5. Associations**

We have performed the same association procedure as in 4FGL, using the Bayesian method and likelihood-ratio method (see Abdollahi et al. 2020). The procedure yielded associations for 341 new DR2 sources and 366 new DR3 sources, representing association fractions of 48% and 41% respectively. These fractions are roughly consistent with that expected from DR1 for these low-TS sources. We provide low-probability (0.1 < P < 0.8) associations for 95 new DR2 or DR3 sources and associations with Planck sources for 44 others. One major change regarding classes is that now we split the generic pulsar class used in previous catalogs into two separate classes, namely young pulsars (PSR) and millisecond pulsars (MSP). We recall that some pulsars with LAT-detected pulsations do not fill the defining significance criterion of this general source catalog and are thus omitted. The full list of LAT-detected pulsars is regularly updated,<fnref n="94" /> and the current population will be discussed in the forthcoming third LAT pulsar catalog. We also point out that associations with pulsar candidates are based on spatial coincidence only and that γ-ray pulsations from these sources are not necessarily expected to be detected soon.

**5.1. Census of New DR2/DR3 Sources Relative to 4FGL-DR1**

For completeness, we recapitulate below the census of associations of new DR2 sources already given in the DR2 document (Ballet et al. 2020), with some updates on blazar classes:

1. two pulsar candidates, PSR J1439−5501 and PSR J1904−11;

2. one globular cluster (GLC), NGC 362;

3. one high-mass binary (HMB), PSR B1259−63; although detected by the LAT during its 2010 periastron passage (Abdo et al. 2011), this source is included in an FGL catalog for the first time thanks to the gain in significance from its 2017 passage (Johnson et al. 2018);

4. three SNRs, 3C 397, SNR G001.4−00.1, and SNR G003.7−00.2; they are considered as reliable associations (as opposed to SPPs, see below) because their radio sizes are quite small (4′−14′) and their positions well match those of the 4FGL sources;

5. two star-forming regions, NGC 346, which is the brightest star-forming region in the SMC, and Sh 2−152;

6. two galaxies, IC 678 and NGC 5380;

7. one starburst galaxy, Arp 299;

8. two radio galaxies, NGC 3078 and NGC 4261;

9. 281 blazars, including 176 blazars of unknown type (BCUs), 66 BL Lacertae objects (BLLs), and 39 FSRQs;<fnref n="95" />

10. 16 SPPs (SPP stands for “SNR, Pulsar or PWN” and refers to sources of unknown nature but overlapping with

known SNRs or PWNe and thus candidate members of these classes);

11. 30 sources of unknown nature (UNK, $|b| < 10°$ sources solely associated with the likelihood-ratio method from large radio and X-ray surveys).

The pulsar PSR J1757−2421, formerly associated with 4FGL J1756.6−2352 in DR1, is now associated with 4FGL J1757.9−2419. The source 4FGL J1756.6−2352 is now a SPP.

Below is a census of the associations for the new DR3 sources:

1. three pulsars with detected pulsations, including two young pulsars, PSR J0922+0638 and PSR J1731−4744 (Smith et al. 2019), and one MSP, PSR J1455−3330;

2. one MSP candidate, PSR J1957+2516;

3. five GLCs, NGC 1851, NGC 5286, NGC 6205, NGC 6712, and M 54;

4. two PWNe, PWN G327.1−1.1 and PWN G54.1+0.3;

5. one SNR, the Kepler SNR; the DR3 source 4FGL J1730.4−2131 was associated by the Likelihood Ratio method with NVSS J173036−212910; this NVSS source is actually part of the Kepler SNR (extended for the Very Large Array at 3′ diameter), so we report the SNR as the counterpart (see also Acero et al. 2022);

6. three novae, V1369 Cen (Cheung et al. 2016), V906 Car (Jean et al. 2018; Aydi et al. 2020), and YZ Ret (MGAB-V207; Li et al. 2020); the position of the DR3 source (4FGL J0358.5−5432) is closer to that of a blazar candidate, PMN J0358−5434, than to that of YZ Ret, but the light curve clearly points to the nova association;

7. one high-mass X-ray binary, 1RXS J172006.1−311702 (Esposito et al. 2014); the nature of this source is debated as it was classified as a cataclysmic variable in Fortin et al. (2018);

8. one radio galaxy, LEDA 55267, classified as a Fanaroff–Riley type 0 (Baldi et al. 2019) with a LAT detection reported in Paliya (2021), and one AGN, NGC 6454;

9. 306 blazar candidates including 219 BCUs, 51 BLLs and 36 FSRQs;

10. 22 SPPs;

11. 20 sources of unknown nature.

**5.2. Association/Classification Changes in 4FGL-DR1 Sources**

1. A specific class (GC) has been created for the Galactic center (4FGL J1745.6−2859), replacing its former SPP class.

2. The increase of localization systematic uncertainties for sources lying at $|b| < 5°$ (Section 3.2) has led to the association of five extra BCUs, five SPPs, and the suppression of the association of 4FGL J1804.9−3001 with the GLC NGC 6528. Four sources classified as unknown in DR1 are now classified as BCUs because they gained Bayesian-based associations. Two other formerly unknown sources are now associated with a pulsar and a binary star (Kleinmann’s star).

3. Two associations with pulsar candidates for 4FGL J1618.7−4633 and 4FGL J1720.8−1937 have been discarded on the realization that their estimated γ-ray luminosities would exceed their spin-down powers by factors greater than 30.

<footnote n="94">⁹⁴ See https://confluence.slac.stanford.edu/display/GLAMCOG/Public+List+of+LAT-Detected+Gamma-Ray+Pulsars.</footnote>

<footnote n="95">⁹⁵ Eight BCUs have been classified into seven BLLs and one FSRQ since the DR2 release.</footnote>
4. Pulsations have been detected for 22 more DR1 sources, including 10 PSR and 12 MSP (e.g., Clark et al. 2021; Tabassum et al. 2021).

5. A total of 29 associations with MSP from the West Virginia University list<fnref n="96" /> have been added (three replacing associations with BCUs and one with a SPP), along with one association with a young pulsar candidate.

6. The pulsar PSR J1909−3744 was mistakenly associated with 4FGL J1912.2−3636. The position of this LAT-detected pulsar was used as a seed, but the resulting 4FGL source ended up too far away from the seed to make the association plausible. The 4FGL source is now unassociated.

7. The association of the MSP binary recently discovered by Wang et al. (2020) with 4FGL J0935.3+0901 has been implemented. Its class is binary (BIN).

8. Following Bogdanov & Halpern (2015), 4FGL J1544.5−1126 is a candidate transitional MSP binary (1RXS J154439.4−112820). Its class, BCU in DR1, has been changed to “low-mass X-ray binary” (LMB). Another candidate transitional MSP binary, CXOU J110926.4−650224, was found (Coti Zelati et al. 2019, 2021) in a systematic search for X-ray counterparts to unassociated sources in the preliminary LAT 8 yr point-source list (FL8Y). While the position of the X-ray source lies well within the 95% error ellipse of the FL8Y source, it is outside of the ellipse of the corresponding DR3 source 4FGL J1110.3−6501. In consequence,<fnref n="97" /> we have not retained this association in our tables.

9. Another transitional MSP binary, 1SXPS J042749.2−670434, displaying simultaneous optical, X-ray, and γ-ray eclipses has been identified by Kennedy et al. (2020) in 4FGL J0427.8−6704 (previously unassociated).

10. Gamma² Velorum (4FGL J0809.5−4714; Martí-Devesa et al. 2020) has been classified as a binary system.

11. We have implemented the results of Corbet et al. (2019), who discovered that 4FGL J1405.1−6119 is a HMB.

12. The source 4FGL J1832.9−0913, associated with the TeV source HESS J1832−093, has been classified as a HMB (see Martí-Devesa & Reimer 2020).

13. Swihart et al. (2021) found a candidate redback MSP binary in a follow-up observation of 4FGL J0940.3−7610. Other redbacks have been discovered in a similar fashion in 4FGL J0407.7−5702 (Miller et al. 2020) and 4FGL J0540.0−7552 (Strader et al. 2021). These associations have been classified as LMBs.

14. Similarly, Li et al. (2021) found a black widow binary at the position of 4FGL J0336.0+7502, which has then also been classified as a LMB.

15. Three periodic optical and X-ray sources were recently discovered positionally consistent with DR1 sources. These redback MSP candidates have been classified as BIN. They are 4FGL J0212.1+5321 (Li et al. 2016), 4FGL J0523.3−2527 (Strader et al. 2014), and 4FGL J0838.7−2827 (Halpern et al. 2017).

16. The tentative association of 4FGL J0647.7−4418 with the HMB RX J0648.0−4418 called out in DR1 has been replaced by association with the BCU SUMSS

J064744−441946 following the multiwavelength investigation of Marti et al. (2020).

17. Following Järvelä et al. (2020), we have reclassified TXS 2116−077 (4FGL J2118.8−0723) as a Seyfert galaxy instead of a NLSY1.

18. The latest version of the Radio Fundamental Catalog<fnref n="98" /> has enabled the association of six previously unassociated sources and two SPPs with blazar candidates.

19. Recent follow-up observations of 4FGL blazars (Peña-Herazo et al. 2017, 2019, 2020, 2021a, 2021b, 2021c; Desai et al. 2019; Rajagopal et al. 2021) have enabled the classification of 240 former BCUs, 2 AGNs, and 2 UNKs into 214 BLLs and 30 FSRQs.

20. Two FSRQs (PKS 0736−770, TXS 1530−131) have been reclassified as BCUs and one other (Two Micron All Sky Survey, 2MASS J02212698+2514338) as a BLL.

Note: A typo has been found in § 6.2 of the 4FGL paper concerning the name of the star cluster associated with the extended H II region encompassing 4FGL J1115.1−6118. The star cluster is NGC 3603 (and not NGC 4603).

Table 5 lists the class tallies of the whole DR3 catalog.

### 5.3. Spectral Properties

In this section, the spectral properties for the main source classes in the whole 4FGL-DR3 sample are presented. We use these properties to shed light on the nature of the unassociated sources in the next section. Although some classes, most notably pulsars, show significant levels of spectral curvature, the PL photon index, Γ, is a convenient parameter enabling comparison across classes. Figure 14 compares the photon-index distributions for different established classes. Pulsars show remarkably consistent soft spectra, along with GLCs. MSP have somewhat harder spectra on average compared to PSR and GLCs. The observed difference between the MSP and GLC distributions seems to challenge the conventional wisdom that MSP dominate the GLC’s γ-ray emission (e.g., Harding et al. 2005; Venter et al. 2009), and is possibly due to different ambient stellar radiation fields comptonized by the MSP electrons (Bednarek & Sitarek 2007). The spectra of SNRs and PWNe are notably harder than those of pulsars, in agreement with the observation that several SNRs and almost all PWNe are detected in the TeV domain. The distribution for blazars spans a broad range of photon index, the distributions for FSRQs and BLLs exhibiting moderate overlap (see Ajello et al. 2020). Blazars constitute the only abundant class showing a hard-spectrum component.

Turning to the spectral curvature, while some blazars show significant curvature, as do most pulsars, the two classes occupy different branches in the LP_TS\_curv versus TS plane (Figure 15 top). These branches merge at the low-TS end as curvature cannot be assessed with high significance in that case. Despite the fact that PSR and MSP show slightly different spectra as can be seen from their photon-index distributions, both fall on the same branch in Figure 15 top.

## 6. Unassociated Sources

The variation of the photon-index distribution of unassociated sources with Galactic latitude indicates that different populations coexist (see Figure 21 of the 4FGL paper), with the

<footnote n="96">⁹⁶ http://astro.phys.wvu.edu/GalacticMSPs/GalacticMSPs.txt</footnote>

<footnote n="97">⁹⁷ Because of the targeted search in the X-rays, applying the standard association procedure would be inappropriate in this case.</footnote>

<footnote n="98">⁹⁸ rfc_2021a available at http://astrogeo.org/rfc/.</footnote>
<table><caption>**Table 5**
LAT 4FGL-DR3 Source Classes</caption><thead><tr><th rowspan="2">Description</th><th colspan="2">Identified</th><th colspan="2">Associated</th></tr><tr><th>Designator</th><th>Number</th><th>Designator</th><th>Number</th></tr></thead><tbody><tr><td>Galactic center</td><td>GC</td><td>1</td><td>...</td><td>...</td></tr><tr><td>Young pulsars, identified by pulsations</td><td>PSR</td><td>135</td><td>...</td><td>...</td></tr><tr><td>Young pulsars, no pulsations seen in LAT yet</td><td>...</td><td>...</td><td>psr</td><td>2</td></tr><tr><td>Millisecond pulsars, identified by pulsations</td><td>MSP</td><td>120</td><td>...</td><td>...</td></tr><tr><td>Millisecond pulsars, no pulsations seen in LAT yet</td><td>...</td><td>...</td><td>msp</td><td>35</td></tr><tr><td>Pulsar wind nebula</td><td>PWN</td><td>11</td><td>pwn</td><td>8</td></tr><tr><td>Supernova remnant</td><td>SNR</td><td>24</td><td>snr</td><td>19</td></tr><tr><td>Supernova remnant / Pulsar wind nebula</td><td>SPP</td><td>0</td><td>spp</td><td>114</td></tr><tr><td>Globular cluster</td><td>GLC</td><td>0</td><td>glc</td><td>35</td></tr><tr><td>Star-forming region</td><td>SFR</td><td>3</td><td>sfr</td><td>2</td></tr><tr><td>High-mass binary</td><td>HMB</td><td>8</td><td>hmb</td><td>3</td></tr><tr><td>Low-mass binary</td><td>LMB</td><td>2</td><td>lmb</td><td>6</td></tr><tr><td>Binary</td><td>BIN</td><td>1</td><td>bin</td><td>6</td></tr><tr><td>Nova</td><td>NOV</td><td>4</td><td>nov</td><td>0</td></tr><tr><td>BL Lac type of blazar</td><td>BLL</td><td>22</td><td>bll</td><td>1435</td></tr><tr><td>FSRQ type of blazar</td><td>FSRQ</td><td>44</td><td>fsrq</td><td>750</td></tr><tr><td>Radio galaxy</td><td>RDG</td><td>6</td><td>rdg</td><td>39</td></tr><tr><td>Nonblazar active galaxy</td><td>AGN</td><td>1</td><td>agn</td><td>8</td></tr><tr><td>Steep spectrum radio quasar</td><td>SSRQ</td><td>0</td><td>ssrq</td><td>2</td></tr><tr><td>Compact steep spectrum radio source</td><td>CSS</td><td>0</td><td>css</td><td>5</td></tr><tr><td>Blazar candidate of uncertain type</td><td>BCU</td><td>1</td><td>bcu</td><td>1491</td></tr><tr><td>Narrow-line Seyfert 1</td><td>NLSY1</td><td>4</td><td>nlsy1</td><td>4</td></tr><tr><td>Seyfert galaxy</td><td>SEY</td><td>0</td><td>sey</td><td>2</td></tr><tr><td>Starburst galaxy</td><td>SBG</td><td>0</td><td>sbg</td><td>8</td></tr><tr><td>Normal galaxy (or part)</td><td>GAL</td><td>2</td><td>gal</td><td>4</td></tr><tr><td>Unknown</td><td>UNK</td><td>0</td><td>unk</td><td>134</td></tr><tr><td>Total</td><td>...</td><td>389</td><td>...</td><td>4112</td></tr><tr><td>Unassociated</td><td>...</td><td>...</td><td>...</td><td>2157</td></tr></tbody></table>

Note. The designation “spp” indicates potential association with SNR or PWN. “Unknown” are $|b| < 10^{\circ}$ sources solely associated with the likelihood-ratio method from large radio and X-ray surveys. Designations shown in capital letters are firm identifications; lower-case letters indicate associations.

spectra becoming softer on average for lower Galactic latitudes. This observation justifies treating separately the cases of high- and low-latitude unassociated sources. We set a somewhat arbitrary limit at $|b| = 10^{\circ}$, with 1122 and 1018 unassociated point sources lying at higher and lower latitudes respectively.

## 6.1. High-latitude Unassociated Sources

Drawing on the class-wise discriminating power of the photon-index distributions illustrated in Figure 14, the distribution for high-latitude unassociated sources is compared to that of BCUs, which are the most relevant blazar set given the low TS, in Figure 16 top. The good agreement between the two distributions supports the conclusion that most high-latitude unassociated sources are indeed blazars.<fnref n="99" /> This

conclusion is slightly qualified when considering the loci of the unassociated sources in the TS versus LP_TScurv plane (Figure 15 bottom). A significant fraction (≈70%) of the higher-TS (TS > 100) unassociated sources fall onto the pulsar branch, strongly suggesting that they are MSP. This idea is also supported by considering the evolution of the photon-index and variability-index distributions with TS (Figure 17). While for the lower-TS bins, the distributions quite closely match those of BCUs, clear differences are observed for the higher-TS bin (100 < TS < 500), where the distributions are compatible with those expected for MSP. In conclusion, their photon indices, significance of spectral curvature, and lack of variability make many of the high-latitude, 100 < TS < 500 unassociated sources good MSP candidates.

## 6.2. Low-latitude Unassociated Sources

Since some AGN catalogs suffer from significant incompleteness close to the plane due to extinction, some background blazars detected through the plane could elude association and contribute to the unassociated sample. Based on the tallies observed at high latitude and the difference in flux limits due to the brighter diffuse emission background close to the plane, the number of $|b| < 10^{\circ}$ blazars is estimated to be 340 ± 20, while 399 blazars are actually reported as associations in DR3. Comparing the low-latitude photon-index distribution with that obtained for higher latitudes (normalized so that the total yield of $\Gamma < 2.4$ sources match), an excess of soft sources ($\Gamma \approx 2.5$) is found in the low-latitude sample relative to the other. From the difference of the spectral-index distributions, this excess, which was less pronounced in 3FGL and is very likely due to a contamination from a nonblazar population, amounts to 75 ± 4 sources. In conclusion, an excess of soft-spectra sources stands out in the low-latitude blazar sample. It can be attributed to a contamination from another population, likely of Galactic origin. After subtracting this excess from the number of reported blazars, the remaining number of detected blazars agrees reasonably well with that extrapolated from the high-latitude population (324 versus 340 ± 20). Few blazars with missed associations are thus expected to contribute to the low-latitude unassociated sample.

Unassociated sources around the Galactic plane (referred to as GUs) account for a growing fraction of the detected sources, from 12% in DR1 to 24% in the new sample. The features presented here were already present in DR1 but have become more salient with the larger population.

The GU latitude distribution, shown in Figure 18, shows a narrow Gaussian-like component (≈210 sources) centered on the plane, and two broad shoulders (≈760 sources) extending out to about $|b| = 10^{\circ}$ (these components will be referred to as “spike” and “shoulders” respectively). The spike is narrower ($\sigma_{\sin(b)} = 0.005$) than the young-pulsar distribution ($\sigma_{\sin(b)} = 0.012$ for sources with $\sin(b) < 0.03$). The loci of the unassociated sources close to the plane exhibit a notable clustering. Some high-density regions (“hot spots”) are readily seen in Figure 19, while other areas are mostly devoid of unassociated sources. The clustering aspect also manifests itself in the large fraction (30%) of GUs with Flag 5 (i.e., close to a brighter neighbor). High-density regions have been determined using a HEALPix tesselation with $N_{side} = 8$ produced in Galactic coordinates. The average density of GUs per pixel in the three pixel rows spanning the Galactic Plane is about 7.

<footnote n="99">⁹⁹ No other source class shows a similar distribution; see Figure 14.</footnote>
<figure><img src="image_16.png" /><figcaption>Figure 17. Variability-index (left) and photon-index (right) distributions of high-latitude unassociated sources and blazars for different bins in Test_Statistic.</figcaption></figure>

<figure><img src="image_17.png" /><figcaption>Figure 18. Galactic-latitude distribution of unassociated sources around the Galactic plane compared to that of pulsars.</figcaption></figure>

most SPP or UNK sources probably belong to the same population as GUs as ascertained by their similar photon-index

distributions (peaking around $\Gamma = 2.5$), they are also flagged if the above conditions are filled. Flag 14 is set for 446 SGUs, 54 SPPs, and 49 UNKs.<fnref n="100" /> SGUs were already significantly present in DR1, but their number in high-density regions has now increased by 50% (see census of Flag 14 in Table 4).

Returning to the spectral properties, the photon-index distribution is essentially constant across the plane, as illustrated in Figure 20 for different (not mutually exclusive) sky regions. In this figure, “clusters” refers to the high-density pixels mentioned above. The unassociated sources preferably lie within the pulsar branch in the LP_TScurv versus TS plane (Figure 21), supporting their spectral resemblance with pulsars. Their distribution of the LP curvature parameter $\beta$ is also consistent with that of pulsars. Despite these similarities, other factors distinguish GUs and pulsars. First, the photon-index distributions do not closely match. The significant clustering

<footnote n="100">¹⁰⁰ We checked that the set of flagged sources depends moderately (typical >85% overlap) on the centering of the HEALPix pixels.</footnote>
<figure><img src="image_18.png" /><figcaption>**Figure 19.** Positions of unassociated sources (crosses) around the Galactic center. The background is a count map of 83–228 MeV *patch* photons simulated over 12 yr with a pixel size of 1°.</figcaption></figure>

<figure><img src="image_19.png" /><figcaption>**Figure 20.** Photon-index distributions of low-latitude unassociated sources in different sky regions.</figcaption></figure>

observed also disfavors the possibility that GUs are mostly pulsars. Moreover, the low-$b$ unassociated sources outnumber the detected PSR by about a factor of 9, making the scenario that so many pulsars have gone undetected so far implausible. The final blow to this possibility comes from their broad latitude distribution, which does not match that of PSR (nor any of the distributions of known γ-ray emitters for that matter). Older pulsars ($>10^6$ yr), which could have drifted farther off the plane over time, exhibit a wider Galactic-latitude distribution, more compatible with that observed. However, these pulsars should be undetectable in γ-rays as they would have crossed the “pair-creation death line” (see Abdo et al. 2013).

Most GUs exhibit a low source-to-background ratio (30% of them having the corresponding Flag 4 set). Whether they originate from mismodeled Galactic diffuse emission becomes a valid question. This emission is not fully understood, and its

modeling suffers many uncertainties. To produce the Pass 8 model,<fnref n="101" /> the contributions of different components, H I, CO, dark neutral medium, and inverse Compton, whose properties are established at larger wavelengths, were fitted to the data in 9 rings of increasing Galactocentric distances. An ad hoc component called the patch encapsulates the contributions of nontemplate components like the Galactic center excess, Loop I, or the Fermi bubbles. Other, yet-to-be-identified large scale structures may be contributing (see Jóhannesson & Porter 2021, for a search performed at high latitude), and a comprehensive analysis of the patch still needs to be carried out. Understanding the patch low-energy behavior would be of particular relevance to this work. The patch covers a limited sky area and was smoothed out not to absorb features with spatial scales less

<footnote n="101">¹⁰¹ https://fermi.gsfc.nasa.gov/ssc/data/analysis/software/aux/4fgl/Galactic_Diffuse_Emission_Model_for_the_4FGL_Catalog_Analysis.pdf</footnote>
<figure><img src="image_20.png" /><figcaption>**Figure 21.** Same as Figure 15 bottom, but for low-latitude unassociated sources.</figcaption></figure>

than 4°. Looking for spatial correlations between the SGUs and features present in the different components, the best match was obtained with the patch component at low energy (see Figure 19). The SGU latitude distribution exhibits a north/south symmetry and a relative extension (Figure 18) that find a good correspondence in the patch at $|b| < 10°$. Comprehending the photon excess underlying the patch will certainly be key in unveiling the nature of the SGUs.

Below we list some possibilities regarding the SGUs and lay out avenues for further exploring them.

1. As mentioned above, SGUs may result from mismodeled Galactic diffuse emission. Allowing more flexibility to the model fitting by leaving free individual diffuse components, like H I and CO, would constrain the robustness of the SGU detections, if the templates used for the model are accurate on scales that are finer than the RoI size, while inaccurate on larger scales.

2. A possibility to be considered is that SGUs trace the large-scale regions associated with fresh injection of cosmic-rays, leading to locally harder γ-ray emission than assumed in the model components. A worthwhile test would be to investigate (e.g., by means of simulations) whether the spurious point-like sources spawned from a possible, not-accounted-for diffuse excess would share common spectral features with the observed SGUs.

3. Some SGUs could be spatially extended if they are related to structures in the diffuse emission. This would introduce some spectral distortion when analyzing them under the assumption that they are point-like, making them softer than they really are.

4. Some correlation with the general structures of interstellar clouds is found by visual inspection, suggesting that the emission ascribed to SGUs could be due to molecular gas unaccounted for on small scales, e.g., because of saturated CO lines. This situation is reminiscent of that encountered in the EGRET era, which led to the discovery of “dark gas” (Grenier et al. 2005) and the suppression of about 30% of unassociated EGRET sources (Casandjian & Grenier 2008) once a modified diffuse emission model was considered. Assessing the effect of saturated CO lines using ¹³CO data would shed light on this possibility.

5. Alternative ways could be pursued to build the patch (diffuse emission that is not accounted for by standard

components of the interstellar emission model) and distinguish it from point sources (e.g., Selig et al. 2015).

6. The alternative to mismodeled diffuse emission, namely the existence of an abundant, entirely new class of sources remains highly speculative. Follow-up multi-wavelength studies (as in Bruzewski et al. 2021), e.g., of regions of high SGU densities or surrounding the brightest SGUs, could help elucidate this issue.

## 7. Conclusions

We built two incremental versions (DR2 after 10 yr and DR3 after 12 yr) of the fourth Fermi-LAT source catalog, the deepest-yet in the GeV energy range. Incremental versions preserve the existence and positions of sources in the previous versions (with a few exceptions, see Section 3.2), and add new sources. We will continue providing incremental versions every two years, until a major ingredient (the event reconstruction, the calibrations, or the underlying diffuse emission model) improves significantly.

We applied a new method to check the quality of the global model of the sky (see Section 3.6). That method indicates that the interstellar emission model is too bright in a number of local molecular clouds. Conversely, it may be too faint close to the Galactic plane where we find many clustered soft unassociated sources (see Section 6.2). Understanding the nature of these sources and what fraction are actually diffuse emission will be an important challenge, which could potentially lead to major advances like the discovery of the dark gas in the EGRET data (Grenier et al. 2005).

The 4FGL-DR3 catalog includes 6658 sources, 1695 of which (25.5%) are found to be significantly variable on one-year timescales, and 3388 of which (50.9%) are fit with curved spectral shapes (see Section 3.4). We mark 317 (4.8%) of the sources as potentially related to imperfections in the model for Galactic diffuse emission; the character c (for “on top of interstellar gas clump”) is appended to their names (except those already marked as e for extended). An additional 1836 (27.6%) are flagged in the catalog for less serious concerns, e.g., for depending sensitively on the details of the analysis or for being close to a brighter source (see Section 3.5).

Of the 6658 sources in the catalog, 389 (5.8%) are considered identified, based on pulsations, correlated variability, or correlated angular sizes with observations at other wavelengths. We find likely lower-energy counterparts for 4112 other sources (61.8%). The remaining 2157 sources (32.4%) are unassociated.

The identified and associated sources in the 4FGL-DR3 catalog include many Galactic and extragalactic source classes (see Section 5). The largest Galactic source class continues to be pulsars, with 135 young, 120 millisecond γ-ray pulsars, and 37 associations to non-LAT pulsars (mostly MSP). Other Galactic source classes have continued to grow; 35 GLCs, 43 SNRs, and 19 PWNe are now associated with LAT sources. Blazars remain the largest class of extragalactic sources, with 2251 identified or associated with BLLs or FSRQ active galaxies and 1492 with unclassified blazars.

The Fermi LAT Collaboration acknowledges generous ongoing support from a number 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
National Aeronautics and Space Administration and the Department of Energy in the United States, the Commissariat à l'Energie Atomique and the Centre National 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 Organization (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 is performed in part under DOE Contract DE-AC02-76SF00515.

This work made extensive use of the ATNF pulsar catalog<sup>102</sup> (Manchester et al. 2005). This research has made use of the NASA/IPAC Extragalactic Database (NED; 2019), which is operated by the Jet Propulsion Laboratory, California Institute of Technology, under contract with the National Aeronautics and Space Administration, and of archival data, software, and online services provided by the ASI Science Data Center (ASDC) operated by the Italian Space Agency. We used the Manitoba SNR catalog (Ferrand & Safi-Harb 2012) to check recently published extended sources. We acknowledge the Einstein@Home project for providing new pulsar associations through the dedicated efforts of the Einstein@Home

volunteers. The Einstein@Home project is supported by the NSF award 1816904.

*Facility:* Fermi-LAT.

*Software:* Gardian (Ackermann et al. 2012), GALPROP<sup>103</sup> (Porter et al. 2017), HEALPix<sup>104</sup> (Górski et al. 2005), Aladin,<sup>105</sup> TOPCAT<sup>106</sup> (Taylor 2005).

# Appendix

## Description of the FITS Version of the 4FGL-DR3 Catalog

The FITS format version of the third release of the 4FGL catalog has seven binary table extensions. The extension LAT_Point_Source_Catalog Extension has all of the information about the sources, and has changed to adapt to the evolutions described in Section 3.4. Its format is described in Table 6. Columns (`Unc_`) `PLEC_Index` and (`Unc_`) `PLEC_Expfactor` have been replaced by (`Unc_`) `PLEC IndexS` and (`Unc_`) `PLEC_ExpfactorS` corresponding to the new pulsar parameterization (Equation (2)). The table has 6659 rows for 6658 sources because the Crab Nebula is described by two entries (the synchrotron component and the inverse Compton component) but counted as only one source. The Crab pulsar is another entry and counted as a separate source. The extended sources are indicated by an “e” appended to their names in the main table.

The extensions GTI, ExtendedSources, ROIs, Components, EnergyBounds, Hist_Start describe auxiliary details and have the same format as in 4FGL.

<footnote n="102">102 http://www.atnf.csiro.au/research/pulsar/psrcat</footnote>

<footnote n="103">103 http://galprop.stanford.edu</footnote>

<footnote n="104">104 http://healpix.jpl.nasa.gov/</footnote>

<footnote n="105">105 http://aladin.u-strasbg.fr/</footnote>

<footnote n="106">106 http://www.star.bristol.ac.uk/~mbt/topcat/</footnote>
<table><caption>Table 6
LAT 4FGL-DR3 FITS Format: LAT_Point_Source_Catalog Extension. Changes since 4FGL Appear in Italics. The File Is Available at https://fermi.gsfc.nasa.gov/ssc/data/access/lat/12yr Catalog/</caption><thead><tr><th>Column</th><th>Format</th><th>Unit</th><th>Description</th></tr></thead><tbody><tr><td>Source_Name</td><td>18A</td><td>...</td><td>Source name 4FGL JHHMM.m+DDMMa<sup>a</sup></td></tr><tr><td>DataRelease</td><td>I</td><td>...</td><td>1 for DR1, 2 for new in DR2, 3 for new or changed in DR3</td></tr><tr><td>RAJ2000</td><td>E</td><td>deg</td><td>Right ascension</td></tr><tr><td>DEJ2000</td><td>E</td><td>deg</td><td>Declination</td></tr><tr><td>GLON</td><td>E</td><td>deg</td><td>Galactic longitude</td></tr><tr><td>GLAT</td><td>E</td><td>deg</td><td>Galactic latitude</td></tr><tr><td>Conf_68_SemiMajor</td><td>E</td><td>deg</td><td>Long radius of error ellipse at 68% confidence<sup>b</sup></td></tr><tr><td>Conf_68_SemiMinor</td><td>E</td><td>deg</td><td>Short radius of error ellipse at 68% confidence<sup>b</sup></td></tr><tr><td>Conf_68_PosAng</td><td>E</td><td>deg</td><td>Position angle of the 68% ellipse<sup>b</sup></td></tr><tr><td>Conf_95_SemiMajor</td><td>E</td><td>deg</td><td>Long radius of error ellipse at 95% confidence</td></tr><tr><td>Conf_95_SemiMinor</td><td>E</td><td>deg</td><td>Short radius of error ellipse at 95% confidence</td></tr><tr><td>Conf_95_PosAng</td><td>E</td><td>deg</td><td>Position angle (eastward) of the long axis from celestial north</td></tr><tr><td>ROI_num</td><td>I</td><td>...</td><td>ROI number (cross-reference to ROIs extension)</td></tr><tr><td>Extended_Source_Name</td><td>18A</td><td>...</td><td>Cross-reference to the ExtendedSources extension</td></tr><tr><td>Signif_Avg</td><td>E</td><td>...</td><td>Source significance in σ units over the 100 MeV to 1 TeV band</td></tr><tr><td>Pivot_Energy</td><td>E</td><td>MeV</td><td>Energy at which error on differential flux is minimal</td></tr><tr><td>Flux1000</td><td>E</td><td>cm<sup>-2</sup> s<sup>-1</sup></td><td>Integral photon flux from 1 to 100 GeV</td></tr><tr><td>Unc_Flux1000</td><td>E</td><td>cm<sup>-2</sup> s<sup>-1</sup></td><td>1σ error on integral photon flux from 1 to 100 GeV</td></tr><tr><td>Energy_Flux100</td><td>E</td><td>erg cm<sup>-2</sup> s<sup>-1</sup></td><td>Energy flux from 100 MeV to 100 GeV obtained by spectral fitting</td></tr><tr><td>Unc_Energy_Flux100</td><td>E</td><td>erg cm<sup>-2</sup> s<sup>-1</sup></td><td>1σ error on energy flux from 100 MeV to 100 GeV</td></tr><tr><td>SpectrumType</td><td>18A</td><td>...</td><td>Spectral type in the global model (PowerLaw, LogParabola, PLSuperExpCutoff)</td></tr><tr><td>PL_Flux_Density</td><td>E</td><td>cm<sup>-2</sup> MeV<sup>-1</sup> s<sup>-1</sup></td><td>Differential flux at Pivot_Energy in PowerLaw fit</td></tr><tr><td>Unc_PL_Flux_Density</td><td>E</td><td>cm<sup>-2</sup> MeV<sup>-1</sup> s<sup>-1</sup></td><td>1σ error on PL_Flux_Density</td></tr><tr><td>PL_Index</td><td>E</td><td>...</td><td>Photon index when fitting with PowerLaw</td></tr><tr><td>Unc_PL_Index</td><td>E</td><td>...</td><td>1σ error on PL_Index</td></tr><tr><td>LP_Flux_Density</td><td>E</td><td>cm<sup>-2</sup> MeV<sup>-1</sup> s<sup>-1</sup></td><td>Differential flux at Pivot_Energy in LogParabola fit</td></tr><tr><td>Unc_LP_Flux_Density</td><td>E</td><td>cm<sup>-2</sup> MeV<sup>-1</sup> s<sup>-1</sup></td><td>1σ error on LP_Flux_Density</td></tr><tr><td>LP_Index</td><td>E</td><td>...</td><td>Photon index at Pivot_Energy (α of Equation (1)) when fitting with LogParabola</td></tr><tr><td>Unc_LP_Index</td><td>E</td><td>...</td><td>1σ error on LP_Index</td></tr><tr><td>LP_beta</td><td>E</td><td>...</td><td>Curvature parameter (β of Equation (1)) when fitting with LogParabola</td></tr><tr><td>Unc_LP_beta</td><td>E</td><td>...</td><td>1σ error on LP_beta</td></tr><tr><td>LP_SigCurv</td><td>E</td><td>...</td><td>Significance (in σ units) of the fit improvement between PowerLaw and LogParabola. A value greater than 4 indicates significant curvature</td></tr><tr><td>LP_EPeak</td><td>E</td><td>MeV</td><td>Peak energy in νF<sub>ν</sub> estimated from the LogParabola model</td></tr><tr><td>Unc_LP_EPeak</td><td>E</td><td>MeV</td><td>1σ error on LP_EPeak</td></tr><tr><td>PLEC_Flux_Density</td><td>E</td><td>cm<sup>-2</sup> MeV<sup>-1</sup> s<sup>-1</sup></td><td>Differential flux at Pivot_Energy in PLSuperExpCutoff fit</td></tr><tr><td>Unc_PLEG_Flux_Density</td><td>E</td><td>cm<sup>-2</sup> MeV<sup>-1</sup> s<sup>-1</sup></td><td>1σ error on PLEG_Flux_Density</td></tr><tr><td>PLEC_IndexS</td><td>E</td><td>...</td><td>Photon index at Pivot_Energy (Γ<sub>s</sub> of Equation (2)) when fitting with PLSuperExpCutoff</td></tr><tr><td>Unc_PLEG_IndexS</td><td>E</td><td>...</td><td>1σ error on PLEG_IndexS</td></tr><tr><td>PLEC_ExpfactorS</td><td>E</td><td>...</td><td>Spectral curvature at Pivot_Energy (d of Equation (2)) when fitting with PLSuperExpCutoff</td></tr><tr><td>Unc_PLEG_ExpfactorS</td><td>E</td><td>...</td><td>1σ error on PLEG_ExpfactorS</td></tr><tr><td>PLEC_Exp_Index</td><td>E</td><td>...</td><td>Exponential index (b of Equation (2)) when fitting with PLSuperExpCutoff</td></tr><tr><td>Unc_PLEG_Exp_Index</td><td>E</td><td>...</td><td>1σ error on PLEG_Exp_Index</td></tr><tr><td>PLEC_SigCurv</td><td>E</td><td>...</td><td>Same as LP_SigCurv for PLSuperExpCutoff model</td></tr><tr><td>PLEC_EPeak</td><td>E</td><td>MeV</td><td>Peak energy in νF<sub>ν</sub> estimated from the PLSuperExpCutoff model</td></tr><tr><td>Unc_PLEG_EPeak</td><td>E</td><td>MeV</td><td>1σ error on PLEG_EPeak</td></tr><tr><td>Npred</td><td>E</td><td>...</td><td>Predicted number of events in the model</td></tr><tr><td>Flux_Band</td><td>8E</td><td>cm<sup>-2</sup> s<sup>-1</sup></td><td>Integral photon flux in each spectral band</td></tr><tr><td>Unc_Flux_Band</td><td>2 × 8E</td><td>cm<sup>-2</sup> s<sup>-1</sup></td><td>1σ lower and upper error on Flux_Band<sup>c</sup></td></tr><tr><td>nuFnu_Band</td><td>8E</td><td>erg cm<sup>-2</sup> s<sup>-1</sup></td><td>Spectral energy distribution over each spectral band</td></tr><tr><td>Sqrt_TS_Band</td><td>8E</td><td>...</td><td>Square root of the Test Statistic in each spectral band</td></tr><tr><td>Variability_Index</td><td>E</td><td>...</td><td>Sum of 2× log(Likelihood) difference between the flux fitted in each time interval and the average flux over the full catalog interval; a value greater than 24.72 over 12 intervals indicates &lt;1% chance of being a steady source</td></tr><tr><td>Frac_Variability</td><td>E</td><td>...</td><td>Fractional variability computed from the fluxes in each year</td></tr><tr><td>Unc_Frac_Variability</td><td>E</td><td>...</td><td>1σ error on fractional variability</td></tr><tr><td>SignifPeak</td><td>E</td><td>...</td><td>Source significance in peak interval in σ units</td></tr><tr><td>Flux_Peak</td><td>E</td><td>cm<sup>-2</sup> s<sup>-1</sup></td><td>Peak integral photon flux from 100 MeV to 100 GeV</td></tr><tr><td>Unc_Flux_Peak</td><td>E</td><td>cm<sup>-2</sup> s<sup>-1</sup></td><td>1σ error on peak integral photon flux</td></tr><tr><td>Time_Peak</td><td>D</td><td>s (MET)</td><td>Time of center of interval in which peak flux was measured</td></tr><tr><td>Peak_Interval</td><td>E</td><td>s</td><td>Length of interval in which peak flux was measured</td></tr></tbody></table>
<table><caption>**Table 6**
(Continued)</caption><thead><tr><th>Column</th><th>Format</th><th>Unit</th><th>Description</th></tr></thead><tbody><tr><td>Flux_History</td><td>12E</td><td>cm<sup>-2</sup> s<sup>-1</sup></td><td>Integral photon flux from 100 MeV to 100 GeV in each year (best fit from likelihood analysis with spectral shape fixed to that obtained over full interval)</td></tr><tr><td>Unc_Flux_History</td><td>2 × 12E</td><td>cm<sup>-2</sup> s<sup>-1</sup></td><td>1σ lower and upper error on integral photon flux in each year<sup>c</sup></td></tr><tr><td>Sqrt_TS_History</td><td>12E</td><td>...</td><td>Square root of the Test Statistic in each year</td></tr><tr><td>ASSOC_4FGL</td><td>18A</td><td>...</td><td>Correspondence to 4FGL source catalog, if any</td></tr><tr><td>ASSOC_FGL</td><td>18A</td><td>...</td><td>Most recent correspondence to FGL source catalogs prior to 4FGL, if any</td></tr><tr><td>ASSOC_FHL</td><td>18A</td><td>...</td><td>Most recent correspondence to previous FHL source catalogs, if any</td></tr><tr><td>ASSOC_GAM1</td><td>18A</td><td>...</td><td>Name of likely corresponding 2AGL source, if any</td></tr><tr><td>ASSOC_GAM2</td><td>18A</td><td>...</td><td>Name of likely corresponding 3EG source, if any</td></tr><tr><td>ASSOC_GAM3</td><td>18A</td><td>...</td><td>Name of likely corresponding EGR source, if any</td></tr><tr><td>TEVCAT_FLAG</td><td>A</td><td>...</td><td>P if positional association with nonextended source in TeVCat<br>E if associated with an extended source in TeVCat, N if no TeV association</td></tr><tr><td>ASSOC_TEV</td><td>24A</td><td>...</td><td>Name of likely corresponding TeV source from TeVCat, if any</td></tr><tr><td>CLASS1</td><td>5A</td><td>...</td><td>Class designation for associated source; see Table 5</td></tr><tr><td>CLASS2</td><td>5A</td><td>...</td><td>Class designation for low-confidence association</td></tr><tr><td>ASSOC1</td><td>28A</td><td>...</td><td>Name of identified or likely associated source</td></tr><tr><td>ASSOC2</td><td>28A</td><td>...</td><td>Name of low-confidence association or of enclosing extended source</td></tr><tr><td>ASSOC_PROB_BAY</td><td>E</td><td>...</td><td>Probability of association according to the Bayesian method<sup>d</sup></td></tr><tr><td>ASSOC_PROB_LR</td><td>E</td><td>...</td><td>Probability of association according to the Likelihood Ratio method<sup>e</sup></td></tr><tr><td>RA_Counterpart</td><td>D</td><td>deg</td><td>R.A. of the counterpart ASSOC1</td></tr><tr><td>DEC_Counterpart</td><td>D</td><td>deg</td><td>decl. of the counterpart ASSOC1</td></tr><tr><td>Unc_Counterpart</td><td>E</td><td>deg</td><td>95% precision of the counterpart localization<sup>f</sup></td></tr><tr><td>Flags</td><td>I</td><td>...</td><td>Source flags (binary coding as in Table 4)<sup>g</sup></td></tr></tbody></table>

**Notes.**

<sup>a</sup> The coordinates are rounded, following the International Astronomical Union convention. The letter at the end can be c (coincident with interstellar clump), e (extended source), i (for Crab Nebula inverse Compton) or s (for Crab Nebula synchrotron).

<sup>b</sup> From the 95% ellipse, assuming a Gaussian distribution.

<sup>c</sup> Separate 1σ errors are computed from the likelihood profile toward lower and larger fluxes. The lower error is set equal to NULL, and the upper error is derived from a Bayesian upper limit if the 1σ interval contains 0 (TS < 1).

<sup>d</sup> NaN in this column when ASSOC1 is defined means that the probability could not be computed, either because the source is extended or because the counterpart is the result of dedicated follow-up.

<sup>e</sup> Probabilities <0.8 are formally set to 0.

<sup>f</sup> For extended counterparts, this reports their extension radius.

<sup>g</sup> Each condition is indicated by one bit among the 16 bits forming Flags. The bit is raised (set to 1) in the dubious case, so that sources without any warning sign have Flags = 0.

**ORCID iDs**

F. Acero https://orcid.org/0000-0002-6606-2816

L. Baldini https://orcid.org/0000-0002-9785-7726

J. Ballet https://orcid.org/0000-0002-8784-2977

D. Bastieri https://orcid.org/0000-0002-6954-8862

R. Bellazzini https://orcid.org/0000-0002-2469-7063

E. Bissaldi https://orcid.org/0000-0001-9935-8106

R. D. Blandford https://orcid.org/0000-0002-1854-5506

R. Bonino https://orcid.org/0000-0002-4264-1215

A. Brill https://orcid.org/0000-0002-6208-5244

R. J. Britto https://orcid.org/0000-0003-3456-2362

P. Bruel https://orcid.org/0000-0002-9032-7941

T. H. Burnett https://orcid.org/0000-0003-0330-9962

R. A. Cameron https://orcid.org/0000-0003-0942-2747

P. A. Caraveo https://orcid.org/0000-0003-2478-8018

D. Castro https://orcid.org/0000-0002-0394-3173

S. Chaty https://orcid.org/0000-0002-5769-8601

C. C. Cheung https://orcid.org/0000-0002-4377-0174

S. Ciprini https://orcid.org/0000-0002-0712-2479

S. Cutini https://orcid.org/0000-0002-1271-2924

F. D'Ammando https://orcid.org/0000-0001-7618-7527

S. De Gaetano https://orcid.org/0000-0002-3358-2559

N. Di Lalla https://orcid.org/0000-0002-7574-1298

L. Di Venere https://orcid.org/0000-0003-0703-824X

A. Domínguez https://orcid.org/0000-0002-3433-4610

V. Fallah Ramazani https://orcid.org/0000-0001-8991-7744

E. C. Ferrara https://orcid.org/0000-0001-7828-7708

A. Fiori https://orcid.org/0000-0003-3174-0688

Y. Fukazawa https://orcid.org/0000-0002-0921-8837

S. Funk https://orcid.org/0000-0002-2012-0080

P. Fusco https://orcid.org/0000-0002-9383-2425

V. Gammaldi https://orcid.org/0000-0003-1826-6117

F. Gargano https://orcid.org/0000-0002-5055-6395

S. Garrappa https://orcid.org/0000-0003-2403-4582

D. Gasparini https://orcid.org/0000-0002-5064-9495

M. Giroletti https://orcid.org/0000-0002-8657-8852

T. Glanzman https://orcid.org/0000-0001-9649-3871

D. Green https://orcid.org/0000-0003-0768-2203

M.-H. Grondin https://orcid.org/0000-0002-8383-251X

L. Guillemot https://orcid.org/0000-0002-9049-8716

S. Guiriec https://orcid.org/0000-0001-5780-8770

A. K. Harding https://orcid.org/0000-0001-6119-859X

E. Hays https://orcid.org/0000-0002-8172-593X

J. W. Hewitt https://orcid.org/0000-0001-5254-2248

X. Hou https://orcid.org/0000-0003-0933-6101

G. Jóhannesson https://orcid.org/0000-0003-1458-7036
M. Kerr https://orcid.org/0000-0002-0893-4073

M. Kuss https://orcid.org/0000-0003-1212-9998

S. Larsson https://orcid.org/0000-0003-0716-107X

L. Latronico https://orcid.org/0000-0002-0984-1856

M. Lemoine-Goumard https://orcid.org/0000-0002-4462-3686

J. Li https://orcid.org/0000-0003-1720-9727

F. Longo https://orcid.org/0000-0003-2501-2270

F. Loparco https://orcid.org/0000-0002-1173-5673

B. Lott https://orcid.org/0000-0003-2186-9242

P. Lubrano https://orcid.org/0000-0003-0221-4806

S. Maldera https://orcid.org/0000-0002-0698-4421

A. Manfreda https://orcid.org/0000-0002-0998-4953

M. N. Mazziotta https://orcid.org/0000-0001-9325-4672

I. Mereu https://orcid.org/0000-0002-0219-4534

N. Mirabal https://orcid.org/0000-0002-7021-5838

W. Mitthumsiri https://orcid.org/0000-0002-3776-072X

T. Mizuno https://orcid.org/0000-0001-7263-0296

M. E. Monzani https://orcid.org/0000-0002-8254-5308

A. Morselli https://orcid.org/0000-0002-7704-9553

I. V. Moskalenko https://orcid.org/0000-0001-6141-458X

E. Nuss https://orcid.org/0000-0002-3528-5777

N. Omodei https://orcid.org/0000-0002-5448-7577

M. Persic https://orcid.org/0000-0003-1853-4900

M. Pesce-Rollins https://orcid.org/0000-0003-1790-8018

V. Petrosian https://orcid.org/0000-0002-2670-8942

R. Pillera https://orcid.org/0000-0003-3808-963X

H. Poon https://orcid.org/0000-0002-6004-4270

T. A. Porter https://orcid.org/0000-0002-2621-4440

G. Principe https://orcid.org/0000-0003-0406-7387

S. Rainò https://orcid.org/0000-0002-9181-0345

R. Rando https://orcid.org/0000-0001-6992-818X

B. Rani https://orcid.org/0000-0001-5711-084X

M. Razzano https://orcid.org/0000-0003-4825-1629

S. Razzaque https://orcid.org/0000-0002-0130-2460

A. Reimer https://orcid.org/0000-0001-8604-7077

O. Reimer https://orcid.org/0000-0001-6953-1385

M. Sánchez-Conde https://orcid.org/0000-0002-3849-9164

P. M. Saz Parkinson https://orcid.org/0000-0001-6566-1246

C. Serini https://orcid.org/0000-0002-9754-6530

C. Sgrò https://orcid.org/0000-0001-5676-6214

D. A. Smith https://orcid.org/0000-0002-7833-0275

D. J. Suson https://orcid.org/0000-0003-2911-2025

D. Tak https://orcid.org/0000-0002-9852-2469

D. J. Thompson https://orcid.org/0000-0001-5217-9135

D. F. Torres https://orcid.org/0000-0002-1522-9065

E. Troja https://orcid.org/0000-0002-1869-7817

J. Valverde https://orcid.org/0000-0002-8090-6528

## References

Abdo, A. A., Ackermann, M., Ajello, M., et al. 2010, A&A, 523, A46

Abdo, A. A., Ackermann, M., Ajello, M., et al. 2011, ApJL, 736, L11

Abdo, A. A., Ajello, M., Allafort, A., et al. 2013, ApJS, 208, 17

Abdollahi, S., Acero, F., Ackermann, M., et al. 2020, ApJS, 247, 33

Abdollahi, S., Ackermann, M., Ajello, M., et al. 2017, ApJ, 846, 34

Acero, F., Lemoine-Goumard, M., & Ballet, J. 2022, A&A, 660, A129

Ackermann, M., Ajello, M., Atwood, W. B., et al. 2012, ApJ, 750, 3

Ajello, M., Angioni, R., Axelsson, M., et al. 2020, ApJ, 892, 105

Ajello, M., Atwood, W. B., Baldini, L., et al. 2017, ApJS, 232, 18

Ambrogi, L., Zanin, R., Casanova, S., et al. 2019, A&A, 623, A86

Araya, M. 2020, MNRAS, 492, 5980

Araya, M., Mitchell, A. M. W., & Parsons, R. D. 2019, MNRAS, 485, 1001

Atwood, W. B., Abdo, A. A., Ackermann, M., et al. 2009, ApJ, 697, 1071

Atwood, W. B., Albert, A., Baldini, L., et al. 2013, arXiv:1303.3514

Aydi, E., Sokolovsky, K. V., Chomiuk, L., et al. 2020, NatAs, 4, 776

Baldi, R. D., Capetti, A., & Giovannini, G. 2019, MNRAS, 482, 2294

Baldini, L., Ballet, J., Bastieri, D., et al. 2021, ApJS, 256, 13

Ballet, J., Burnett, T. H., Digel, S. W., & Lott, B. 2020, arXiv:2005.11208

Bednarek, W., & Sitarek, J. 2007, MNRAS, 377, 920

Bogdanov, S., & Halpern, J. P. 2015, ApJL, 803, L27

Bruel, P. 2021, A&A, 656, A81

Bruel, P., Burnett, T. H., Digel, S. W., et al. 2018, arXiv:1810.11394

Bruzewski, S., Schinzel, F. K., Taylor, G. B., & Petrov, L. 2021, ApJ, 914, 42

Casandjian, J.-M., & Grenier, I. A. 2008, A&A, 489, 849

Cheung, C. C., Jean, P., Shore, S. N., et al. 2016, ApJ, 826, 142

Clark, C. J., Nieder, L., Voisin, G., et al. 2021, MNRAS, 502, 915

Corbet, R. H. D., Chomiuk, L., Coe, M. J., et al. 2019, ApJ, 884, 93

Coti Zelati, F., Hugo, B., Torres, D. F., et al. 2021, A&A, 655, A52

Coti Zelati, F., Papitto, A., de Martino, D., et al. 2019, A&A, 622, A211

Desai, A., Marchesi, S., Rajagopal, M., & Ajello, M. 2019, ApJS, 241, 5

Devin, J., Lemoine-Goumard, M., Grondin, M. H., et al. 2020, A&A, 643, A28

Esposito, P., Israel, G. L., Sidoli, L., et al. 2014, MNRAS, 441, 1126

Ferrand, G., & Safi-Harb, S. 2012, AdSpR, 49, 1313

Fortin, F., Chaty, S., Coleiro, A., Tomsick, J. A., & Nitschelm, C. H. R. 2018, A&A, 618, A150

Górski, K. M., Hivon, E., Banday, A. J., et al. 2005, ApJ, 622, 759

Grenier, I. A., Casandjian, J.-M., & Terrier, R. 2005, Sci, 307, 1292

Grondin, M.-H., Funk, S., Lemoine-Goumard, M., et al. 2011, ApJ, 738, 42

Halpern, J. P., Strader, J., & Li, M. 2017, ApJ, 844, 150

Harding, A. K., Usov, V. V., & Muslimov, A. G. 2005, ApJ, 622, 531

Järvelä, E., Berton, M., Ciroi, S., et al. 2020, A&A, 636, L12

Jean, P., Cheung, C. C., Ojha, R., van Zyl, P., & Angioni, R. 2018, ATel, 11546, 1

Jóhannesson, G., & Porter, T. A. 2021, ApJ, 917, 30

Johnson, T. J., Wood, K. S., Kerr, M., et al. 2018, ApJ, 863, 27

Kennedy, M. R., Breton, R. P., Clark, C. J., et al. 2020, MNRAS, 494, 3912

Kerr, M. 2010, PhD Thesis, University of Washington

Li, K.-L., Jane Yap, Y. X., Hui, C. Y., & Kong, A. K. H. 2021, ApJ, 911, 92

Li, K.-L., Kong, A., Aydi, E., et al. 2020, ATel, 13868, 1

Li, K.-L., Kong, A. K. H., Hou, X., et al. 2016, ApJ, 833, 143

Manchester, R. N., Hobbs, G. B., Teoh, A., & Hobbs, M. 2005, AJ, 129, 1993

Mares, A., Lemoine-Goumard, M., Acero, F., et al. 2021, ApJ, 912, 158

Marti, J., Sanchez-Ayaso, E., Luque-Escamilla, P. L., et al. 2020, MNRAS, 492, 4291

Martí-Devesa, G., & Reimer, O. 2020, A&A, 637, A23

Martí-Devesa, G., Reimer, O., Li, J., & Torres, D. F. 2020, A&A, 635, A141

McConville, W., Ostorero, L., Moderski, R., et al. 2011, ApJ, 738, 148

Miller, J. M., Swihart, S. J., Strader, J., et al. 2020, ApJ, 904, 49

NASA/IPAC Extragalactic Database (NED) 2019, NASA/IPAC Extragalactic Database (NED), IPAC, doi:10.26132/NED1

Paliya, V. S. 2021, ApJL, 918, L39

Peña-Herazo, H. A., Amaya-Almazán, R. A., Massaro, F., et al. 2020, A&A, 643, A103

Peña-Herazo, H. A., Marchesini, E. J., Álvarez Crespo, N., et al. 2017, Ap&SS, 362, 228

Peña-Herazo, H. A., Massaro, F., Chavushyan, V., et al. 2019, Ap&SS, 364, 85

Peña-Herazo, H. A., Massaro, F., Gu, M., et al. 2021b, AJ, 162, 76

Peña-Herazo, H. A., Massaro, F., Gu, M., et al. 2021c, AJ, 161, 196

Peña-Herazo, H. A., Paggi, A., García-Pérez, A., et al. 2021a, AJ, 162, 177

Porter, T. A., Jóhannesson, G., & Moskalenko, I. V. 2017, ApJ, 846, 67

Principe, G., Mitchell, A. M. W., Caroff, S., et al. 2020, A&A, 640, A76

Rajagopal, M., Marchesi, S., Kaur, A., et al. 2021, ApJS, 254, 26

Selig, M., Vacca, V., Oppermann, N., & Enßlin, T. A. 2015, A&A, 581, A126

Smith, D. A., Bruel, P., Cognard, I., et al. 2019, ApJ, 871, 78

Strader, J., Chomiuk, L., Sonbas, E., et al. 2014, ApJL, 788, L27

Strader, J., Swihart, S. J., Urquhart, R., et al. 2021, ApJ, 917, 69

Swihart, S. J., Strader, J., Aydi, E., et al. 2021, ApJ, 909, 185

Tabassum, S., Ransom, S., Ray, P., et al. 2021, 43rd COSPAR Scientific Assembly, 1208

Taylor, M. B. 2005, in ASP Conf. Ser., 347, Astronomical Data Analysis Software and Systems XIV, ed. P. Shopbell, M. Britton, & R. Ebert (San Francisco, CA: ASP), 29

Tibaldo, L., Zanin, R., Faggioli, G., et al. 2018, A&A, 617, A78

Venter, C., De Jager, O. C., & Clapson, A. C. 2009, ApJL, 696, L52

Wang, Z., Xing, Y., Zhang, J., et al. 2020, MNRAS, 493, 4845

Xin, Y., Zeng, H., Liu, S., Fan, Y., & Wei, D. 2019, ApJ, 885, 162