Cosmic Ray Spectrum and Composition from PeV to EeV Using 3 Years of Data From IceTop and IceCube

Cosmic Ray Spectrum and Composition from PeV to EeV Using 3 Years of Data From IceTop and IceCube

III. Physikalisches Institut, RWTH Aachen University, D-52056 Aachen, Germany Department of Physics, University of Adelaide, Adelaide, 5005, Australia Dept. of Physics and Astronomy, University of Alaska Anchorage, 3211 Providence Dr., Anchorage, AK 99508, USA Dept. of Physics, University of Texas at Arlington, 502 Yates St., Science Hall Rm 108, Box 19059, Arlington, TX 76019, USA CTSPS, Clark-Atlanta University, Atlanta, GA 30314, USA School of Physics and Center for Relativistic Astrophysics, Georgia Institute of Technology, Atlanta, GA 30332, USA Dept. of Physics, Southern University, Baton Rouge, LA 70813, USA Dept. of Physics, University of California, Berkeley, CA 94720, USA Lawrence Berkeley National Laboratory, Berkeley, CA 94720, USA Institut für Physik, Humboldt-Universität zu Berlin, D-12489 Berlin, Germany Fakultät für Physik & Astronomie, Ruhr-Universität Bochum, D-44780 Bochum, Germany Université Libre de Bruxelles, Science Faculty CP230, B-1050 Brussels, Belgium Vrije Universiteit Brussel (VUB), Dienst ELEM, B-1050 Brussels, Belgium Dept. of Physics, Massachusetts Institute of Technology, Cambridge, MA 02139, USA Dept. of Physics and Institute for Global Prominent Research, Chiba University, Chiba 263-8522, Japan Dept. of Physics and Astronomy, University of Canterbury, Private Bag 4800, Christchurch, New Zealand Dept. of Physics, University of Maryland, College Park, MD 20742, USA Dept. of Astronomy, Ohio State University, Columbus, OH 43210, USA Dept. of Physics and Center for Cosmology and Astro-Particle Physics, Ohio State University, Columbus, OH 43210, USA Niels Bohr Institute, University of Copenhagen, DK-2100 Copenhagen, Denmark Dept. of Physics, TU Dortmund University, D-44221 Dortmund, Germany Dept. of Physics and Astronomy, Michigan State University, East Lansing, MI 48824, USA Dept. of Physics, University of Alberta, Edmonton, Alberta, Canada T6G 2E1 Erlangen Centre for Astroparticle Physics, Friedrich-Alexander-Universität Erlangen-Nürnberg, D-91058 Erlangen, Germany Physik-department, Technische Universität München, D-85748 Garching, Germany Département de physique nucléaire et corpusculaire, Université de Genève, CH-1211 Genève, Switzerland Dept. of Physics and Astronomy, University of Gent, B-9000 Gent, Belgium Dept. of Physics and Astronomy, University of California, Irvine, CA 92697, USA Dept. of Physics and Astronomy, University of Kansas, Lawrence, KS 66045, USA SNOLAB, 1039 Regional Road 24, Creighton Mine 9, Lively, ON, Canada P3Y 1N2 Department of Physics and Astronomy, UCLA, Los Angeles, CA 90095, USA Dept. of Astronomy, University of Wisconsin, Madison, WI 53706, USA Dept. of Physics and Wisconsin IceCube Particle Astrophysics Center, University of Wisconsin, Madison, WI 53706, USA Institute of Physics, University of Mainz, Staudinger Weg 7, D-55099 Mainz, Germany Department of Physics, Marquette University, Milwaukee, WI, 53201, USA Institut für Kernphysik, Westfälische Wilhelms-Universität Münster, D-48149 Münster, Germany Bartol Research Institute and Dept. of Physics and Astronomy, University of Delaware, Newark, DE 19716, USA Dept. of Physics, Yale University, New Haven, CT 06520, USA Dept. of Physics, University of Oxford, Parks Road, Oxford OX1 3PQ, UK Dept. of Physics, Drexel University, 3141 Chestnut Street, Philadelphia, PA 19104, USA Physics Department, South Dakota School of Mines and Technology, Rapid City, SD 57701, USA Dept. of Physics, University of Wisconsin, River Falls, WI 54022, USA Dept. of Physics and Astronomy, University of Rochester, Rochester, NY 14627, USA Oskar Klein Centre and Dept. of Physics, Stockholm University, SE-10691 Stockholm, Sweden Dept. of Physics and Astronomy, Stony Brook University, Stony Brook, NY 11794-3800, USA Dept. of Physics, Sungkyunkwan University, Suwon 16419, Korea Dept. of Physics and Astronomy, University of Alabama, Tuscaloosa, AL 35487, USA Dept. of Astronomy and Astrophysics, Pennsylvania State University, University Park, PA 16802, USA Dept. of Physics, Pennsylvania State University, University Park, PA 16802, USA Dept. of Physics and Astronomy, Uppsala University, Box 516, S-75120 Uppsala, Sweden Dept. of Physics, University of Wuppertal, D-42119 Wuppertal, Germany DESY, D-15738 Zeuthen, Germany    M. G. Aartsen Dept. of Physics and Astronomy, University of Canterbury, Private Bag 4800, Christchurch, New Zealand    M. Ackermann DESY, D-15738 Zeuthen, Germany    J. Adams Dept. of Physics and Astronomy, University of Canterbury, Private Bag 4800, Christchurch, New Zealand    J. A. Aguilar Université Libre de Bruxelles, Science Faculty CP230, B-1050 Brussels, Belgium    M. Ahlers Niels Bohr Institute, University of Copenhagen, DK-2100 Copenhagen, Denmark    M. Ahrens Oskar Klein Centre and Dept. of Physics, Stockholm University, SE-10691 Stockholm, Sweden    C. Alispach Département de physique nucléaire et corpusculaire, Université de Genève, CH-1211 Genève, Switzerland    K. Andeen Department of Physics, Marquette University, Milwaukee, WI, 53201, USA    T. Anderson Dept. of Physics, Pennsylvania State University, University Park, PA 16802, USA    I. Ansseau Université Libre de Bruxelles, Science Faculty CP230, B-1050 Brussels, Belgium    G. Anton Erlangen Centre for Astroparticle Physics, Friedrich-Alexander-Universität Erlangen-Nürnberg, D-91058 Erlangen, Germany    C. Argüelles Dept. of Physics, Massachusetts Institute of Technology, Cambridge, MA 02139, USA    J. Auffenberg III. Physikalisches Institut, RWTH Aachen University, D-52056 Aachen, Germany    S. Axani Dept. of Physics, Massachusetts Institute of Technology, Cambridge, MA 02139, USA    P. Backes III. Physikalisches Institut, RWTH Aachen University, D-52056 Aachen, Germany    H. Bagherpour Dept. of Physics and Astronomy, University of Canterbury, Private Bag 4800, Christchurch, New Zealand    X. Bai Physics Department, South Dakota School of Mines and Technology, Rapid City, SD 57701, USA    A. Barbano Département de physique nucléaire et corpusculaire, Université de Genève, CH-1211 Genève, Switzerland    S. W. Barwick Dept. of Physics and Astronomy, University of California, Irvine, CA 92697, USA    V. Baum Institute of Physics, University of Mainz, Staudinger Weg 7, D-55099 Mainz, Germany    S. Baur Université Libre de Bruxelles, Science Faculty CP230, B-1050 Brussels, Belgium    R. Bay Dept. of Physics, University of California, Berkeley, CA 94720, USA    J. J. Beatty Dept. of Physics and Center for Cosmology and Astro-Particle Physics, Ohio State University, Columbus, OH 43210, USA Dept. of Astronomy, Ohio State University, Columbus, OH 43210, USA    K.-H. Becker Dept. of Physics, University of Wuppertal, D-42119 Wuppertal, Germany    J. Becker Tjus Fakultät für Physik & Astronomie, Ruhr-Universität Bochum, D-44780 Bochum, Germany    S. BenZvi Dept. of Physics and Astronomy, University of Rochester, Rochester, NY 14627, USA    D. Berley Dept. of Physics, University of Maryland, College Park, MD 20742, USA    E. Bernardini DESY, D-15738 Zeuthen, Germany    D. Z. Besson Dept. of Physics and Astronomy, University of Kansas, Lawrence, KS 66045, USA    G. Binder Lawrence Berkeley National Laboratory, Berkeley, CA 94720, USA Dept. of Physics, University of California, Berkeley, CA 94720, USA    D. Bindig Dept. of Physics, University of Wuppertal, D-42119 Wuppertal, Germany    E. Blaufuss Dept. of Physics, University of Maryland, College Park, MD 20742, USA    S. Blot DESY, D-15738 Zeuthen, Germany    C. Bohm Oskar Klein Centre and Dept. of Physics, Stockholm University, SE-10691 Stockholm, Sweden    M. Börner Dept. of Physics, TU Dortmund University, D-44221 Dortmund, Germany    S. Böser Institute of Physics, University of Mainz, Staudinger Weg 7, D-55099 Mainz, Germany    O. Botner Dept. of Physics and Astronomy, Uppsala University, Box 516, S-75120 Uppsala, Sweden    J. Böttcher III. Physikalisches Institut, RWTH Aachen University, D-52056 Aachen, Germany    E. Bourbeau Niels Bohr Institute, University of Copenhagen, DK-2100 Copenhagen, Denmark    J. Bourbeau Dept. of Physics and Wisconsin IceCube Particle Astrophysics Center, University of Wisconsin, Madison, WI 53706, USA    F. Bradascio DESY, D-15738 Zeuthen, Germany    J. Braun Dept. of Physics and Wisconsin IceCube Particle Astrophysics Center, University of Wisconsin, Madison, WI 53706, USA    H.-P. Bretz DESY, D-15738 Zeuthen, Germany    S. Bron Département de physique nucléaire et corpusculaire, Université de Genève, CH-1211 Genève, Switzerland    J. Brostean-Kaiser DESY, D-15738 Zeuthen, Germany    A. Burgman Dept. of Physics and Astronomy, Uppsala University, Box 516, S-75120 Uppsala, Sweden    J. Buscher III. Physikalisches Institut, RWTH Aachen University, D-52056 Aachen, Germany    R. S. Busse Dept. of Physics and Wisconsin IceCube Particle Astrophysics Center, University of Wisconsin, Madison, WI 53706, USA    T. Carver Département de physique nucléaire et corpusculaire, Université de Genève, CH-1211 Genève, Switzerland    C. Chen School of Physics and Center for Relativistic Astrophysics, Georgia Institute of Technology, Atlanta, GA 30332, USA    E. Cheung Dept. of Physics, University of Maryland, College Park, MD 20742, USA    D. Chirkin Dept. of Physics and Wisconsin IceCube Particle Astrophysics Center, University of Wisconsin, Madison, WI 53706, USA    K. Clark SNOLAB, 1039 Regional Road 24, Creighton Mine 9, Lively, ON, Canada P3Y 1N2    L. Classen Institut für Kernphysik, Westfälische Wilhelms-Universität Münster, D-48149 Münster, Germany    G. H. Collin Dept. of Physics, Massachusetts Institute of Technology, Cambridge, MA 02139, USA    J. M. Conrad Dept. of Physics, Massachusetts Institute of Technology, Cambridge, MA 02139, USA    P. Coppin Vrije Universiteit Brussel (VUB), Dienst ELEM, B-1050 Brussels, Belgium    P. Correa Vrije Universiteit Brussel (VUB), Dienst ELEM, B-1050 Brussels, Belgium    D. F. Cowen Dept. of Physics, Pennsylvania State University, University Park, PA 16802, USA Dept. of Astronomy and Astrophysics, Pennsylvania State University, University Park, PA 16802, USA    R. Cross Dept. of Physics and Astronomy, University of Rochester, Rochester, NY 14627, USA    P. Dave School of Physics and Center for Relativistic Astrophysics, Georgia Institute of Technology, Atlanta, GA 30332, USA    J. P. A. M. de André Dept. of Physics and Astronomy, Michigan State University, East Lansing, MI 48824, USA    C. De Clercq Vrije Universiteit Brussel (VUB), Dienst ELEM, B-1050 Brussels, Belgium    J. J. DeLaunay Dept. of Physics, Pennsylvania State University, University Park, PA 16802, USA    H. Dembinski Bartol Research Institute and Dept. of Physics and Astronomy, University of Delaware, Newark, DE 19716, USA    K. Deoskar Oskar Klein Centre and Dept. of Physics, Stockholm University, SE-10691 Stockholm, Sweden    S. De Ridder Dept. of Physics and Astronomy, University of Gent, B-9000 Gent, Belgium    P. Desiati Dept. of Physics and Wisconsin IceCube Particle Astrophysics Center, University of Wisconsin, Madison, WI 53706, USA    K. D. de Vries Vrije Universiteit Brussel (VUB), Dienst ELEM, B-1050 Brussels, Belgium    G. de Wasseige Vrije Universiteit Brussel (VUB), Dienst ELEM, B-1050 Brussels, Belgium    M. de With Institut für Physik, Humboldt-Universität zu Berlin, D-12489 Berlin, Germany    T. DeYoung Dept. of Physics and Astronomy, Michigan State University, East Lansing, MI 48824, USA    A. Diaz Dept. of Physics, Massachusetts Institute of Technology, Cambridge, MA 02139, USA    J. C. Díaz-Vélez Dept. of Physics and Wisconsin IceCube Particle Astrophysics Center, University of Wisconsin, Madison, WI 53706, USA    H. Dujmovic Dept. of Physics, Sungkyunkwan University, Suwon 16419, Korea    M. Dunkman Dept. of Physics, Pennsylvania State University, University Park, PA 16802, USA    E. Dvorak Physics Department, South Dakota School of Mines and Technology, Rapid City, SD 57701, USA    B. Eberhardt Dept. of Physics and Wisconsin IceCube Particle Astrophysics Center, University of Wisconsin, Madison, WI 53706, USA    T. Ehrhardt Institute of Physics, University of Mainz, Staudinger Weg 7, D-55099 Mainz, Germany    P. Eller Dept. of Physics, Pennsylvania State University, University Park, PA 16802, USA    P. A. Evenson Bartol Research Institute and Dept. of Physics and Astronomy, University of Delaware, Newark, DE 19716, USA    S. Fahey Dept. of Physics and Wisconsin IceCube Particle Astrophysics Center, University of Wisconsin, Madison, WI 53706, USA    A. R. Fazely Dept. of Physics, Southern University, Baton Rouge, LA 70813, USA    J. Felde Dept. of Physics, University of Maryland, College Park, MD 20742, USA    T. Feusels Dept. of Physics and Astronomy, University of Gent, B-9000 Gent, Belgium    K. Filimonov Dept. of Physics, University of California, Berkeley, CA 94720, USA    C. Finley Oskar Klein Centre and Dept. of Physics, Stockholm University, SE-10691 Stockholm, Sweden    A. Franckowiak DESY, D-15738 Zeuthen, Germany    E. Friedman Dept. of Physics, University of Maryland, College Park, MD 20742, USA    A. Fritz Institute of Physics, University of Mainz, Staudinger Weg 7, D-55099 Mainz, Germany    T. K. Gaisser Bartol Research Institute and Dept. of Physics and Astronomy, University of Delaware, Newark, DE 19716, USA    J. Gallagher Dept. of Astronomy, University of Wisconsin, Madison, WI 53706, USA    E. Ganster III. Physikalisches Institut, RWTH Aachen University, D-52056 Aachen, Germany    S. Garrappa DESY, D-15738 Zeuthen, Germany    L. Gerhardt Lawrence Berkeley National Laboratory, Berkeley, CA 94720, USA    K. Ghorbani Dept. of Physics and Wisconsin IceCube Particle Astrophysics Center, University of Wisconsin, Madison, WI 53706, USA    T. Glauch Physik-department, Technische Universität München, D-85748 Garching, Germany    T. Glüsenkamp Erlangen Centre for Astroparticle Physics, Friedrich-Alexander-Universität Erlangen-Nürnberg, D-91058 Erlangen, Germany    A. Goldschmidt Lawrence Berkeley National Laboratory, Berkeley, CA 94720, USA    J. G. Gonzalez Bartol Research Institute and Dept. of Physics and Astronomy, University of Delaware, Newark, DE 19716, USA    D. Grant Dept. of Physics and Astronomy, Michigan State University, East Lansing, MI 48824, USA    Z. Griffith Dept. of Physics and Wisconsin IceCube Particle Astrophysics Center, University of Wisconsin, Madison, WI 53706, USA    M. Günder III. Physikalisches Institut, RWTH Aachen University, D-52056 Aachen, Germany    M. Gündüz Fakultät für Physik & Astronomie, Ruhr-Universität Bochum, D-44780 Bochum, Germany    C. Haack III. Physikalisches Institut, RWTH Aachen University, D-52056 Aachen, Germany    A. Hallgren Dept. of Physics and Astronomy, Uppsala University, Box 516, S-75120 Uppsala, Sweden    L. Halve III. Physikalisches Institut, RWTH Aachen University, D-52056 Aachen, Germany    F. Halzen Dept. of Physics and Wisconsin IceCube Particle Astrophysics Center, University of Wisconsin, Madison, WI 53706, USA    K. Hanson Dept. of Physics and Wisconsin IceCube Particle Astrophysics Center, University of Wisconsin, Madison, WI 53706, USA    D. Hebecker Institut für Physik, Humboldt-Universität zu Berlin, D-12489 Berlin, Germany    D. Heereman Université Libre de Bruxelles, Science Faculty CP230, B-1050 Brussels, Belgium    P. Heix III. Physikalisches Institut, RWTH Aachen University, D-52056 Aachen, Germany    K. Helbing Dept. of Physics, University of Wuppertal, D-42119 Wuppertal, Germany    R. Hellauer Dept. of Physics, University of Maryland, College Park, MD 20742, USA    F. Henningsen Physik-department, Technische Universität München, D-85748 Garching, Germany    S. Hickford Dept. of Physics, University of Wuppertal, D-42119 Wuppertal, Germany    J. Hignight Dept. of Physics and Astronomy, Michigan State University, East Lansing, MI 48824, USA    G. C. Hill Department of Physics, University of Adelaide, Adelaide, 5005, Australia    K. D. Hoffman Dept. of Physics, University of Maryland, College Park, MD 20742, USA    R. Hoffmann Dept. of Physics, University of Wuppertal, D-42119 Wuppertal, Germany    T. Hoinka Dept. of Physics, TU Dortmund University, D-44221 Dortmund, Germany    B. Hokanson-Fasig Dept. of Physics and Wisconsin IceCube Particle Astrophysics Center, University of Wisconsin, Madison, WI 53706, USA    K. Hoshina Dept. of Physics and Wisconsin IceCube Particle Astrophysics Center, University of Wisconsin, Madison, WI 53706, USA    F. Huang Dept. of Physics, Pennsylvania State University, University Park, PA 16802, USA    M. Huber Physik-department, Technische Universität München, D-85748 Garching, Germany    K. Hultqvist Oskar Klein Centre and Dept. of Physics, Stockholm University, SE-10691 Stockholm, Sweden    M. Hünnefeld Dept. of Physics, TU Dortmund University, D-44221 Dortmund, Germany    R. Hussain Dept. of Physics and Wisconsin IceCube Particle Astrophysics Center, University of Wisconsin, Madison, WI 53706, USA    S. In Dept. of Physics, Sungkyunkwan University, Suwon 16419, Korea    N. Iovine Université Libre de Bruxelles, Science Faculty CP230, B-1050 Brussels, Belgium    A. Ishihara Dept. of Physics and Institute for Global Prominent Research, Chiba University, Chiba 263-8522, Japan    E. Jacobi DESY, D-15738 Zeuthen, Germany    G. S. Japaridze CTSPS, Clark-Atlanta University, Atlanta, GA 30314, USA    M. Jeong Dept. of Physics, Sungkyunkwan University, Suwon 16419, Korea    K. Jero Dept. of Physics and Wisconsin IceCube Particle Astrophysics Center, University of Wisconsin, Madison, WI 53706, USA    B. J. P. Jones Dept. of Physics, University of Texas at Arlington, 502 Yates St., Science Hall Rm 108, Box 19059, Arlington, TX 76019, USA    F. Jonske III. Physikalisches Institut, RWTH Aachen University, D-52056 Aachen, Germany    R. Joppe III. Physikalisches Institut, RWTH Aachen University, D-52056 Aachen, Germany    W. Kang Dept. of Physics, Sungkyunkwan University, Suwon 16419, Korea    A. Kappes Institut für Kernphysik, Westfälische Wilhelms-Universität Münster, D-48149 Münster, Germany    D. Kappesser Institute of Physics, University of Mainz, Staudinger Weg 7, D-55099 Mainz, Germany    T. Karg DESY, D-15738 Zeuthen, Germany    M. Karl Physik-department, Technische Universität München, D-85748 Garching, Germany    A. Karle Dept. of Physics and Wisconsin IceCube Particle Astrophysics Center, University of Wisconsin, Madison, WI 53706, USA    U. Katz Erlangen Centre for Astroparticle Physics, Friedrich-Alexander-Universität Erlangen-Nürnberg, D-91058 Erlangen, Germany    M. Kauer Dept. of Physics and Wisconsin IceCube Particle Astrophysics Center, University of Wisconsin, Madison, WI 53706, USA    J. L. Kelley Dept. of Physics and Wisconsin IceCube Particle Astrophysics Center, University of Wisconsin, Madison, WI 53706, USA    A. Kheirandish Dept. of Physics and Wisconsin IceCube Particle Astrophysics Center, University of Wisconsin, Madison, WI 53706, USA    J. Kim Dept. of Physics, Sungkyunkwan University, Suwon 16419, Korea    T. Kintscher DESY, D-15738 Zeuthen, Germany    J. Kiryluk Dept. of Physics and Astronomy, Stony Brook University, Stony Brook, NY 11794-3800, USA    T. Kittler Erlangen Centre for Astroparticle Physics, Friedrich-Alexander-Universität Erlangen-Nürnberg, D-91058 Erlangen, Germany    S. R. Klein Lawrence Berkeley National Laboratory, Berkeley, CA 94720, USA Dept. of Physics, University of California, Berkeley, CA 94720, USA    R. Koirala Bartol Research Institute and Dept. of Physics and Astronomy, University of Delaware, Newark, DE 19716, USA    H. Kolanoski Institut für Physik, Humboldt-Universität zu Berlin, D-12489 Berlin, Germany    L. Köpke Institute of Physics, University of Mainz, Staudinger Weg 7, D-55099 Mainz, Germany    C. Kopper Dept. of Physics and Astronomy, Michigan State University, East Lansing, MI 48824, USA    S. Kopper Dept. of Physics and Astronomy, University of Alabama, Tuscaloosa, AL 35487, USA    D. J. Koskinen Niels Bohr Institute, University of Copenhagen, DK-2100 Copenhagen, Denmark    M. Kowalski Institut für Physik, Humboldt-Universität zu Berlin, D-12489 Berlin, Germany DESY, D-15738 Zeuthen, Germany    K. Krings Physik-department, Technische Universität München, D-85748 Garching, Germany    G. Krückl Institute of Physics, University of Mainz, Staudinger Weg 7, D-55099 Mainz, Germany    N. Kulacz Dept. of Physics, University of Alberta, Edmonton, Alberta, Canada T6G 2E1    S. Kunwar DESY, D-15738 Zeuthen, Germany    N. Kurahashi Dept. of Physics, Drexel University, 3141 Chestnut Street, Philadelphia, PA 19104, USA    A. Kyriacou Department of Physics, University of Adelaide, Adelaide, 5005, Australia    M. Labare Dept. of Physics and Astronomy, University of Gent, B-9000 Gent, Belgium    J. L. Lanfranchi Dept. of Physics, Pennsylvania State University, University Park, PA 16802, USA    M. J. Larson Dept. of Physics, University of Maryland, College Park, MD 20742, USA    F. Lauber Dept. of Physics, University of Wuppertal, D-42119 Wuppertal, Germany    J. P. Lazar Dept. of Physics and Wisconsin IceCube Particle Astrophysics Center, University of Wisconsin, Madison, WI 53706, USA    K. Leonard Dept. of Physics and Wisconsin IceCube Particle Astrophysics Center, University of Wisconsin, Madison, WI 53706, USA    M. Leuermann III. Physikalisches Institut, RWTH Aachen University, D-52056 Aachen, Germany    Q. R. Liu Dept. of Physics and Wisconsin IceCube Particle Astrophysics Center, University of Wisconsin, Madison, WI 53706, USA    E. Lohfink Institute of Physics, University of Mainz, Staudinger Weg 7, D-55099 Mainz, Germany    C. J. Lozano Mariscal Institut für Kernphysik, Westfälische Wilhelms-Universität Münster, D-48149 Münster, Germany    L. Lu Dept. of Physics and Institute for Global Prominent Research, Chiba University, Chiba 263-8522, Japan    F. Lucarelli Département de physique nucléaire et corpusculaire, Université de Genève, CH-1211 Genève, Switzerland    J. Lünemann Vrije Universiteit Brussel (VUB), Dienst ELEM, B-1050 Brussels, Belgium    W. Luszczak Dept. of Physics and Wisconsin IceCube Particle Astrophysics Center, University of Wisconsin, Madison, WI 53706, USA    J. Madsen Dept. of Physics, University of Wisconsin, River Falls, WI 54022, USA    G. Maggi Vrije Universiteit Brussel (VUB), Dienst ELEM, B-1050 Brussels, Belgium    K. B. M. Mahn Dept. of Physics and Astronomy, Michigan State University, East Lansing, MI 48824, USA    Y. Makino Dept. of Physics and Institute for Global Prominent Research, Chiba University, Chiba 263-8522, Japan    P. Mallik III. Physikalisches Institut, RWTH Aachen University, D-52056 Aachen, Germany    K. Mallot Dept. of Physics and Wisconsin IceCube Particle Astrophysics Center, University of Wisconsin, Madison, WI 53706, USA    S. Mancina Dept. of Physics and Wisconsin IceCube Particle Astrophysics Center, University of Wisconsin, Madison, WI 53706, USA    I. C. Mariş Université Libre de Bruxelles, Science Faculty CP230, B-1050 Brussels, Belgium    R. Maruyama Dept. of Physics, Yale University, New Haven, CT 06520, USA    K. Mase Dept. of Physics and Institute for Global Prominent Research, Chiba University, Chiba 263-8522, Japan    R. Maunu Dept. of Physics, University of Maryland, College Park, MD 20742, USA    K. Meagher Dept. of Physics and Wisconsin IceCube Particle Astrophysics Center, University of Wisconsin, Madison, WI 53706, USA    M. Medici Niels Bohr Institute, University of Copenhagen, DK-2100 Copenhagen, Denmark    A. Medina Dept. of Physics and Center for Cosmology and Astro-Particle Physics, Ohio State University, Columbus, OH 43210, USA    M. Meier Dept. of Physics, TU Dortmund University, D-44221 Dortmund, Germany    S. Meighen-Berger Physik-department, Technische Universität München, D-85748 Garching, Germany    T. Menne Dept. of Physics, TU Dortmund University, D-44221 Dortmund, Germany    G. Merino Dept. of Physics and Wisconsin IceCube Particle Astrophysics Center, University of Wisconsin, Madison, WI 53706, USA    T. Meures Université Libre de Bruxelles, Science Faculty CP230, B-1050 Brussels, Belgium    S. Miarecki Lawrence Berkeley National Laboratory, Berkeley, CA 94720, USA Dept. of Physics, University of California, Berkeley, CA 94720, USA    J. Micallef Dept. of Physics and Astronomy, Michigan State University, East Lansing, MI 48824, USA    G. Momenté Institute of Physics, University of Mainz, Staudinger Weg 7, D-55099 Mainz, Germany    T. Montaruli Département de physique nucléaire et corpusculaire, Université de Genève, CH-1211 Genève, Switzerland    R. W. Moore Dept. of Physics, University of Alberta, Edmonton, Alberta, Canada T6G 2E1    R. Morse Dept. of Physics and Wisconsin IceCube Particle Astrophysics Center, University of Wisconsin, Madison, WI 53706, USA    M. Moulai Dept. of Physics, Massachusetts Institute of Technology, Cambridge, MA 02139, USA    P. Muth III. Physikalisches Institut, RWTH Aachen University, D-52056 Aachen, Germany    R. Nagai Dept. of Physics and Institute for Global Prominent Research, Chiba University, Chiba 263-8522, Japan    R. Nahnhauer DESY, D-15738 Zeuthen, Germany    P. Nakarmi Dept. of Physics and Astronomy, University of Alabama, Tuscaloosa, AL 35487, USA    U. Naumann Dept. of Physics, University of Wuppertal, D-42119 Wuppertal, Germany    G. Neer Dept. of Physics and Astronomy, Michigan State University, East Lansing, MI 48824, USA    H. Niederhausen Physik-department, Technische Universität München, D-85748 Garching, Germany    S. C. Nowicki Dept. of Physics, University of Alberta, Edmonton, Alberta, Canada T6G 2E1    D. R. Nygren Lawrence Berkeley National Laboratory, Berkeley, CA 94720, USA    A. Obertacke Pollmann Dept. of Physics, University of Wuppertal, D-42119 Wuppertal, Germany    A. Olivas Dept. of Physics, University of Maryland, College Park, MD 20742, USA    A. O’Murchadha Université Libre de Bruxelles, Science Faculty CP230, B-1050 Brussels, Belgium    E. O’Sullivan Oskar Klein Centre and Dept. of Physics, Stockholm University, SE-10691 Stockholm, Sweden    T. Palczewski Lawrence Berkeley National Laboratory, Berkeley, CA 94720, USA Dept. of Physics, University of California, Berkeley, CA 94720, USA    H. Pandya Bartol Research Institute and Dept. of Physics and Astronomy, University of Delaware, Newark, DE 19716, USA    D. V. Pankova Dept. of Physics, Pennsylvania State University, University Park, PA 16802, USA    N. Park Dept. of Physics and Wisconsin IceCube Particle Astrophysics Center, University of Wisconsin, Madison, WI 53706, USA    P. Peiffer Institute of Physics, University of Mainz, Staudinger Weg 7, D-55099 Mainz, Germany    C. Pérez de los Heros Dept. of Physics and Astronomy, Uppsala University, Box 516, S-75120 Uppsala, Sweden    S. Philippen III. Physikalisches Institut, RWTH Aachen University, D-52056 Aachen, Germany    D. Pieloth Dept. of Physics, TU Dortmund University, D-44221 Dortmund, Germany    E. Pinat Université Libre de Bruxelles, Science Faculty CP230, B-1050 Brussels, Belgium    A. Pizzuto Dept. of Physics and Wisconsin IceCube Particle Astrophysics Center, University of Wisconsin, Madison, WI 53706, USA    M. Plum Department of Physics, Marquette University, Milwaukee, WI, 53201, USA    A. Porcelli Dept. of Physics and Astronomy, University of Gent, B-9000 Gent, Belgium    P. B. Price Dept. of Physics, University of California, Berkeley, CA 94720, USA    G. T. Przybylski Lawrence Berkeley National Laboratory, Berkeley, CA 94720, USA    C. Raab Université Libre de Bruxelles, Science Faculty CP230, B-1050 Brussels, Belgium    A. Raissi Dept. of Physics and Astronomy, University of Canterbury, Private Bag 4800, Christchurch, New Zealand    M. Rameez Niels Bohr Institute, University of Copenhagen, DK-2100 Copenhagen, Denmark    L. Rauch DESY, D-15738 Zeuthen, Germany    K. Rawlins Dept. of Physics and Astronomy, University of Alaska Anchorage, 3211 Providence Dr., Anchorage, AK 99508, USA    I. C. Rea Physik-department, Technische Universität München, D-85748 Garching, Germany    R. Reimann III. Physikalisches Institut, RWTH Aachen University, D-52056 Aachen, Germany    B. Relethford Dept. of Physics, Drexel University, 3141 Chestnut Street, Philadelphia, PA 19104, USA    G. Renzi Université Libre de Bruxelles, Science Faculty CP230, B-1050 Brussels, Belgium    E. Resconi Physik-department, Technische Universität München, D-85748 Garching, Germany    W. Rhode Dept. of Physics, TU Dortmund University, D-44221 Dortmund, Germany    M. Richman Dept. of Physics, Drexel University, 3141 Chestnut Street, Philadelphia, PA 19104, USA    S. Robertson Lawrence Berkeley National Laboratory, Berkeley, CA 94720, USA    M. Rongen III. Physikalisches Institut, RWTH Aachen University, D-52056 Aachen, Germany    C. Rott Dept. of Physics, Sungkyunkwan University, Suwon 16419, Korea    T. Ruhe Dept. of Physics, TU Dortmund University, D-44221 Dortmund, Germany    D. Ryckbosch Dept. of Physics and Astronomy, University of Gent, B-9000 Gent, Belgium    D. Rysewyk Dept. of Physics and Astronomy, Michigan State University, East Lansing, MI 48824, USA    I. Safa Dept. of Physics and Wisconsin IceCube Particle Astrophysics Center, University of Wisconsin, Madison, WI 53706, USA    S. E. Sanchez Herrera Dept. of Physics, University of Alberta, Edmonton, Alberta, Canada T6G 2E1    A. Sandrock Dept. of Physics, TU Dortmund University, D-44221 Dortmund, Germany    J. Sandroos Institute of Physics, University of Mainz, Staudinger Weg 7, D-55099 Mainz, Germany    M. Santander Dept. of Physics and Astronomy, University of Alabama, Tuscaloosa, AL 35487, USA    S. Sarkar Dept. of Physics, University of Oxford, Parks Road, Oxford OX1 3PQ, UK    S. Sarkar Dept. of Physics, University of Alberta, Edmonton, Alberta, Canada T6G 2E1    K. Satalecka DESY, D-15738 Zeuthen, Germany    M. Schaufel III. Physikalisches Institut, RWTH Aachen University, D-52056 Aachen, Germany    P. Schlunder Dept. of Physics, TU Dortmund University, D-44221 Dortmund, Germany    T. Schmidt Dept. of Physics, University of Maryland, College Park, MD 20742, USA    A. Schneider Dept. of Physics and Wisconsin IceCube Particle Astrophysics Center, University of Wisconsin, Madison, WI 53706, USA    J. Schneider Erlangen Centre for Astroparticle Physics, Friedrich-Alexander-Universität Erlangen-Nürnberg, D-91058 Erlangen, Germany    L. Schumacher III. Physikalisches Institut, RWTH Aachen University, D-52056 Aachen, Germany    S. Sclafani Dept. of Physics, Drexel University, 3141 Chestnut Street, Philadelphia, PA 19104, USA    D. Seckel Bartol Research Institute and Dept. of Physics and Astronomy, University of Delaware, Newark, DE 19716, USA    S. Seunarine Dept. of Physics, University of Wisconsin, River Falls, WI 54022, USA    S. Shefali III. Physikalisches Institut, RWTH Aachen University, D-52056 Aachen, Germany    M. Silva Dept. of Physics and Wisconsin IceCube Particle Astrophysics Center, University of Wisconsin, Madison, WI 53706, USA    R. Snihur Dept. of Physics and Wisconsin IceCube Particle Astrophysics Center, University of Wisconsin, Madison, WI 53706, USA    J. Soedingrekso Dept. of Physics, TU Dortmund University, D-44221 Dortmund, Germany    D. Soldin Bartol Research Institute and Dept. of Physics and Astronomy, University of Delaware, Newark, DE 19716, USA    M. Song Dept. of Physics, University of Maryland, College Park, MD 20742, USA    G. M. Spiczak Dept. of Physics, University of Wisconsin, River Falls, WI 54022, USA    C. Spiering DESY, D-15738 Zeuthen, Germany    J. Stachurska DESY, D-15738 Zeuthen, Germany    M. Stamatikos Dept. of Physics and Center for Cosmology and Astro-Particle Physics, Ohio State University, Columbus, OH 43210, USA    T. Stanev Bartol Research Institute and Dept. of Physics and Astronomy, University of Delaware, Newark, DE 19716, USA    A. Stasik DESY, D-15738 Zeuthen, Germany    R. Stein DESY, D-15738 Zeuthen, Germany    J. Stettner III. Physikalisches Institut, RWTH Aachen University, D-52056 Aachen, Germany    A. Steuer Institute of Physics, University of Mainz, Staudinger Weg 7, D-55099 Mainz, Germany    T. Stezelberger Lawrence Berkeley National Laboratory, Berkeley, CA 94720, USA    R. G. Stokstad Lawrence Berkeley National Laboratory, Berkeley, CA 94720, USA    A. Stößl Dept. of Physics and Institute for Global Prominent Research, Chiba University, Chiba 263-8522, Japan    N. L. Strotjohann DESY, D-15738 Zeuthen, Germany    T. Stürwald III. Physikalisches Institut, RWTH Aachen University, D-52056 Aachen, Germany    T. Stuttard Niels Bohr Institute, University of Copenhagen, DK-2100 Copenhagen, Denmark    G. W. Sullivan Dept. of Physics, University of Maryland, College Park, MD 20742, USA    M. Sutherland Dept. of Physics and Center for Cosmology and Astro-Particle Physics, Ohio State University, Columbus, OH 43210, USA    I. Taboada School of Physics and Center for Relativistic Astrophysics, Georgia Institute of Technology, Atlanta, GA 30332, USA    F. Tenholt Fakultät für Physik & Astronomie, Ruhr-Universität Bochum, D-44780 Bochum, Germany    S. Ter-Antonyan Dept. of Physics, Southern University, Baton Rouge, LA 70813, USA    A. Terliuk DESY, D-15738 Zeuthen, Germany    S. Tilav Bartol Research Institute and Dept. of Physics and Astronomy, University of Delaware, Newark, DE 19716, USA    L. Tomankova Fakultät für Physik & Astronomie, Ruhr-Universität Bochum, D-44780 Bochum, Germany    C. Tönnis Dept. of Physics, Sungkyunkwan University, Suwon 16419, Korea    S. Toscano Vrije Universiteit Brussel (VUB), Dienst ELEM, B-1050 Brussels, Belgium    D. Tosi Dept. of Physics and Wisconsin IceCube Particle Astrophysics Center, University of Wisconsin, Madison, WI 53706, USA    M. Tselengidou Erlangen Centre for Astroparticle Physics, Friedrich-Alexander-Universität Erlangen-Nürnberg, D-91058 Erlangen, Germany    C. F. Tung School of Physics and Center for Relativistic Astrophysics, Georgia Institute of Technology, Atlanta, GA 30332, USA    A. Turcati Physik-department, Technische Universität München, D-85748 Garching, Germany    R. Turcotte III. Physikalisches Institut, RWTH Aachen University, D-52056 Aachen, Germany    C. F. Turley Dept. of Physics, Pennsylvania State University, University Park, PA 16802, USA    B. Ty Dept. of Physics and Wisconsin IceCube Particle Astrophysics Center, University of Wisconsin, Madison, WI 53706, USA    E. Unger Dept. of Physics and Astronomy, Uppsala University, Box 516, S-75120 Uppsala, Sweden    M. A. Unland Elorrieta Institut für Kernphysik, Westfälische Wilhelms-Universität Münster, D-48149 Münster, Germany    M. Usner DESY, D-15738 Zeuthen, Germany    J. Vandenbroucke Dept. of Physics and Wisconsin IceCube Particle Astrophysics Center, University of Wisconsin, Madison, WI 53706, USA    W. Van Driessche Dept. of Physics and Astronomy, University of Gent, B-9000 Gent, Belgium    D. van Eijk Dept. of Physics and Wisconsin IceCube Particle Astrophysics Center, University of Wisconsin, Madison, WI 53706, USA    N. van Eijndhoven Vrije Universiteit Brussel (VUB), Dienst ELEM, B-1050 Brussels, Belgium    S. Vanheule Dept. of Physics and Astronomy, University of Gent, B-9000 Gent, Belgium    J. van Santen DESY, D-15738 Zeuthen, Germany    M. Vraeghe Dept. of Physics and Astronomy, University of Gent, B-9000 Gent, Belgium    C. Walck Oskar Klein Centre and Dept. of Physics, Stockholm University, SE-10691 Stockholm, Sweden    A. Wallace Department of Physics, University of Adelaide, Adelaide, 5005, Australia    M. Wallraff III. Physikalisches Institut, RWTH Aachen University, D-52056 Aachen, Germany    N. Wandkowsky Dept. of Physics and Wisconsin IceCube Particle Astrophysics Center, University of Wisconsin, Madison, WI 53706, USA    T. B. Watson Dept. of Physics, University of Texas at Arlington, 502 Yates St., Science Hall Rm 108, Box 19059, Arlington, TX 76019, USA    C. Weaver Dept. of Physics, University of Alberta, Edmonton, Alberta, Canada T6G 2E1    M. J. Weiss Dept. of Physics, Pennsylvania State University, University Park, PA 16802, USA    J. Weldert Institute of Physics, University of Mainz, Staudinger Weg 7, D-55099 Mainz, Germany    C. Wendt Dept. of Physics and Wisconsin IceCube Particle Astrophysics Center, University of Wisconsin, Madison, WI 53706, USA    J. Werthebach Dept. of Physics and Wisconsin IceCube Particle Astrophysics Center, University of Wisconsin, Madison, WI 53706, USA    S. Westerhoff Dept. of Physics and Wisconsin IceCube Particle Astrophysics Center, University of Wisconsin, Madison, WI 53706, USA    B. J. Whelan Department of Physics, University of Adelaide, Adelaide, 5005, Australia    N. Whitehorn Department of Physics and Astronomy, UCLA, Los Angeles, CA 90095, USA    K. Wiebe Institute of Physics, University of Mainz, Staudinger Weg 7, D-55099 Mainz, Germany    C. H. Wiebusch III. Physikalisches Institut, RWTH Aachen University, D-52056 Aachen, Germany    L. Wille Dept. of Physics and Wisconsin IceCube Particle Astrophysics Center, University of Wisconsin, Madison, WI 53706, USA    D. R. Williams Dept. of Physics and Astronomy, University of Alabama, Tuscaloosa, AL 35487, USA    L. Wills Dept. of Physics, Drexel University, 3141 Chestnut Street, Philadelphia, PA 19104, USA    M. Wolf Physik-department, Technische Universität München, D-85748 Garching, Germany    J. Wood Dept. of Physics and Wisconsin IceCube Particle Astrophysics Center, University of Wisconsin, Madison, WI 53706, USA    T. R. Wood Dept. of Physics, University of Alberta, Edmonton, Alberta, Canada T6G 2E1    K. Woschnagg Dept. of Physics, University of California, Berkeley, CA 94720, USA    G. Wrede Erlangen Centre for Astroparticle Physics, Friedrich-Alexander-Universität Erlangen-Nürnberg, D-91058 Erlangen, Germany    D. L. Xu Dept. of Physics and Wisconsin IceCube Particle Astrophysics Center, University of Wisconsin, Madison, WI 53706, USA    X. W. Xu Dept. of Physics, Southern University, Baton Rouge, LA 70813, USA    Y. Xu Dept. of Physics and Astronomy, Stony Brook University, Stony Brook, NY 11794-3800, USA    J. P. Yanez Dept. of Physics, University of Alberta, Edmonton, Alberta, Canada T6G 2E1    G. Yodh Dept. of Physics and Astronomy, University of California, Irvine, CA 92697, USA    S. Yoshida Dept. of Physics and Institute for Global Prominent Research, Chiba University, Chiba 263-8522, Japan    T. Yuan Dept. of Physics and Wisconsin IceCube Particle Astrophysics Center, University of Wisconsin, Madison, WI 53706, USA    M. Zöcklein III. Physikalisches Institut, RWTH Aachen University, D-52056 Aachen, Germany analysis@icecube.wisc.edu
July 1, 2019:  VERSION 2.65
Abstract

We report on measurements of the all-particle cosmic ray energy spectrum and composition in the PeV to EeV energy range using three years of data from the IceCube Neutrino Observatory. The IceTop detector measures cosmic ray induced air showers on the surface of the ice, from which the energy spectrum of cosmic rays is determined by making additional assumptions about the mass composition. A separate measurement is performed when IceTop data are analyzed in coincidence with the high-energy muon energy loss information from the deep in-ice IceCube detector. In this measurement, both the spectrum and the mass composition of the primary cosmic rays are simultaneously reconstructed using a neural network trained on observables from both detectors. The performance and relative advantages of these two distinct analyses are discussed, including the systematic uncertainties and the dependence on the hadronic interaction models, and both all-particle spectra as well as individual spectra for elemental groups are presented.

pacs:
thanks: also at Università di Padova, I-35131 Padova, Italythanks: Earthquake Research Institute, University of Tokyo, Bunkyo, Tokyo 113-0032, Japan

IceCube Collaboration

Measurements of the cosmic ray energy spectrum and mass composition in the PeV to EeV energy range provide key information about the origin and propagation of cosmic rays in what is commonly considered to be a transition region from galactic to extragalactic cosmic rays. In the present paper we report on composition and energy spectrum measurements using three years of data from the IceCube Neutrino Observatory.

The IceCube Neutrino Observatory is a versatile particle detector located at the geographic South Pole, with both surface and deeply-buried components. The latter, or InIce, component (described in detail in Aartsen et al. (2017)) consists of 5160 Digital Optical Modules (DOMs) Abbasi et al. (2008a) deployed with 17 m spacing on 86 strings in a 125 m triangular grid formation, at depths from 1450 m to 2450 m below the surface in transparent ice. Each DOM contains a 10 inch Hamamatsu photomultiplier tube (PMT) and electronics for signal processing and readout Abbasi et al. (2009).

The surface component, IceTop (described in detail in Abbasi et al. (2013a)), is an array of pairs of tanks filled with water that has frozen and containing two DOMs each, operating at different PMT gains for increased dynamic range; the pairs of tanks are called stations, and each station is located above a string of the InIce detector (see Figure 1).

As Cherenkov detectors, both components are sensitive to the charged particles in the extensive air showers (EAS) produced by cosmic rays from the Southern Hemisphere sky, and both are used to extract the cosmic ray energy spectrum in different ways. In this paper, we will discuss two methods for extracting this information: the first uses IceTop only and will be referred to as the IceTop-alone analysis, while the other uses both the IceTop and InIce detectors in tandem, and will be referred to as the coincident analysis. The coincident analysis also provides a new measurement of the composition of cosmic rays.

The IceTop detector can be used alone to measure the core position, direction, and size of the air showers at the surface. These observables are utilized to extract a cosmic ray energy spectrum: in particular the shower size, discussed in Section I.4, is strongly correlated with the energy of the incident cosmic ray primary. It is important to note that as the shower size depends slightly on the mass of the cosmic ray primary, the composition must be assumed in order to extract the energy spectrum when using IceTop-alone. Since IceTop collects a large number of showers, rare high-energy events are collected in sufficient numbers to extend the analysis into the EeV range. This analysis was performed most recently in Aartsen et al. (2013a), in which one year of 73-station data was used to measure a spectrum from a few PeV to EeV.

The InIce detector measures the energy loss of the high-energy muons in the deep ice, which is strongly dependent upon the mass of the incident cosmic ray primary, as discussed in Section I.5. Thus, with the IceTop and InIce detectors working in tandem, the high-energy muon component of the air showers are measured in coincidence with the electromagnetic component; thus both the energy spectrum and mass composition are measured without making assumptions about one to determine the other. However, this method requires coincident events between the two components of IceCube: due to the long lever arm between the two arrays, only zenith angles of approximately 0-30^{\circ} on the sky pass the coincident selection criteria, which then yield fewer events from the same data sample as the IceTop-alone analysis. Thus, the coincident analysis cannot reach as high an energy as IceTop-alone. A coincident analysis like this was performed in Abbasi et al. (2013b), in which one month of data from the half-completed 40-station, 40-string detector was used to measure a spectrum and average logarithmic mass (\langle\ln A\rangle) from 1 to 30 PeV; that analysis was extended with improved reconstruction techniques in Feusels et al. (2013); Feusels (2014) using one year of data from the nearly complete 73-station, 79-string array, achieving better resolution and reaching to 1 EeV.

In this paper, we will present for the first time the coincident analysis in detail. In addition, the one-year analyses of IceTop-alone Aartsen et al. (2013a) and coincident Feusels et al. (2013); Feusels (2014) data are extended to three years, and we report these improved and updated results.

I Data, Simulation, and Reconstruction

I.1 The 3-year Data Set

Figure 1: A top view of the IceTop surface array. Colors indicate the construction periods for the strings and tanks. This work will focus on IceTop-73 (IT-73) and IceCube-79 (IC-79), which includes all detectors except those from the last construction year colored orange. Strings 79 and 80 and stations 79, 80, and 81 are excluded from this work, even though they lie inside the border of IT-73 shown in black. Feusels (2014)

The analyses described here use three years of data, from June 1, 2010, through May 2, 2013. The IceTop-alone and the coincident analyses use the same dataset, and thus have both a total livetime of 977.6 days (with a negligible uncertainty of less than half a percent). The IceCube Neutrino Observatory was not running in its complete configuration (81 stations, 86 strings) until May 14, 2011. Therefore, the first year of data included here was taken in the incomplete 73-station, 79-string configuration (IT-73/IC-79). The simulation used in these analyses was also produced using the IT-73/IC-79 detector configuration. Thus, in order to handle the two following years of data in the same analysis, all information from the final 7 deployed strings and 8 stations has been removed from data processing. The configuration used is shown in Figure 1.

I.2 Simulation

Both of the analyses presented here use the same set of Monte Carlo simulations of cosmic ray events to establish relationships between detector observables and cosmic ray energy and mass. (This simulated dataset was also used in Aartsen et al. (2013a).) 30000 air showers of four primary types (protons, helium, oxygen, and iron) were simulated using the CORSIKA air shower generator Heck et al. (1998) with an E^{-1} spectrum between \log_{10}(E/\mathrm{GeV}) = 5.0 and 8.0. Additionally, thinned CORSIKA showers Heck and Pierog (); Hillas (1997) were generated at higher energies: 12000 showers of each type, again with an E^{-1} spectrum between \log_{10}(E/\mathrm{GeV}) = 7.0 and 9.5. High-energy muons are not thinned in the algorithm. The range of energy overlap between \log_{10}(E/\mathrm{GeV}) = 7.0 and 8.0 allows for verification of the un-thinning algorithm Abbasi et al. (2013a).

Two hadronic interaction models were used: FLUKA Battistoni et al. (2007) below 80 GeV, and SYBILL 2.1 Ahn et al. (2009) above 80 GeV. EGS4 Nelson et al. (1985) was used to model the electromagnetic interactions. Other high energy hadronic interaction models are used for systematic studies, as discussed in Section IV.3. The zenith angle from the primary particles is sampled from a \cos(\theta)\sin(\theta) distribution between 0{}^{\circ} and 40{}^{\circ}, while the azimuth angle is drawn from a uniform distribution over the whole 2\pi azimuth range. To be able to study atmospheric changes and apply corrections where needed, a reference atmosphere was chosen based on the MSIS-90-E parametrization Hedin (1991) of the South Pole Atmosphere on July 1, 1997, which has a ground pressure of 692.2 g/cm{}^{2} at the South Pole altitude (CORSIKA atmosphere 12, Heck and Pierog ()). To make more efficient use of the CORSIKA showers available, each shower is copied, or resampled, 100 times, and thrown at random locations within a circle of radius R centered on the center of the IceTop array. The resampling radius is the largest possible for the shower to trigger the array Abbasi et al. (2013a).

Next, the particles generated by CORSIKA are propagated into the detectors. The individual responses of the IceTop tanks are simulated using a detailed Geant4 Agostinelli et al. (2003); Allison et al. (2006) model which takes into account the individual IceTop tank properties, including snow and air above the tanks, and the detector electronics. The resulting signal is converted into units of photoelectrons using constants unique to each tank Abbasi et al. (2013a). In this way the same calibration procedure (described in Section I.3) can be performed on both simulated events and experimental data.

For the coincident analysis, the high-energy muons in the CORSIKA air showers must also be propagated through the Antarctic ice and through the deep InIce detector. Muons with energy above 273 GeV111273 GeV is the energy at which 0.1% of the muons are expected to reach the top of the InIce detector and could create a detectable amount of Cherenkov light within the array Feusels (2014) are propagated Chirkin and Rhode (2004) through the Antarctic ice to the bottom of the InIce array. Propagating the Cherenkov photons from the muons through the South Pole ice to the DOMs by directly tracking each one is computationally prohibitive. Therefore, light profiles for GeV emitters (including both pure Cherenkov emission and more diffuse emission from cascade light sources) are tabulated in a software package Lundberg (2007), which includes a model of the full structure of the ice properties Aartsen et al. (2013b). The expected number of photoelectrons and their arrival times at each DOM are retrieved from the tables for each muon and for electromagnetic and hadronic cascades. Simulated noise hits are then added. Finally, the simulated photoelectrons are fed into a simulation of the InIce readout electronics and the detector trigger. The simulated DOM signals then follow the same processing chain as the experimentally measured DOM signals.

I.3 Pulse cleaning and Calibration

Only events which pass the cosmic ray filters Abbasi et al. (2013a) and which contain at least 6 hard local coincidence (HLC222HLC hits occur when both tanks in one station are hit within 1 \mus) DOMs within 6 \mus are processed further. In a first signal cleaning of selected events, signals from all DOMs which were determined to be unreliable at the time of data taking are removed. The remaining signals are calibrated using the procedure described in Aartsen et al. (2017), which returns the number of photoelectrons. At this point the IceTop and InIce data are split.

The IceTop data is cleaned and calibrated as described in detail in Abbasi et al. (2013a). At this stage, only HLC hits are preserved. All signals within one large trigger window are then split into clusters of hits that are likely related to one air shower. Next, an additional noise cleaning procedure is applied to the calibrated HLC signals. This procedure begins with a cluster of 3 hit stations and adds more tanks to the event using a simple distance-per-time requirement Feusels (2014). This extra cleaning procedure improves the removal of noise hits, of tanks with time-fluctuations, and of muon-like signals at large lateral distances from the shower core. These, in turn, improve the IceTop reconstruction procedures described in Section I.4. The specific tank response is then taken into account using a calibration procedure which transforms the signal into that expected from vertical (equivalent) muons (VEMs).

For the coincident analysis, the InIce hits causally connected to the event seen in IceTop are selected based on an allowed time difference window between IceTop and InIce hits. Therefore, before the IceTop and InIce events can be connected, InIce noise hits that affect the trigger time must first be removed in a procedure based on distance-per-time requirements (similar to that in IceTop). Then, InIce triggers matching the IceTop triggers and pulses are selected, and the pulses not connected to the selected trigger are removed. A final noise removal procedure in the InIce detector uses the reconstruction of the track by IceTop whereby only hits connected in time to the track, and in a cylinder around it, are preserved. These steps effectively remove random coincidences between the two detectors, and events that are very close in time.

Finally, after the above event cleaning and calibration procedure, events are required to have hits in at least 5 IceTop (IT-73) stations and, in the coincident analysis, hits in at least 8 InIce (IC-79) DOMs.

I.4 IceTop Reconstruction

Cleaned data from IceTop tanks are processed by a reconstruction software package called Laputop, which has been described in detail in Abbasi et al. (2013a). For each event, Laputop finds the best-fit shower core position (x_{c}, y_{c}, z_{c}) and direction (\theta, \phi), as well as two parameters describing the shape of the lateral distribution function (LDF) of deposited charge (S_{125}, \beta). The functional form of the charge LDF is a double logarithmic parabola:

S=S_{\mathrm{125}}\cdot\left(\frac{r}{125\,m}\right)^{-\beta-\kappa\log_{10}% \left(\frac{r}{125\,m}\right)},

where \beta is a measure of the steepness of the LDF and S_{125} is the signal expectation at a perpendicular reference distance of 125 meters from the shower axis Klepser (2008). S_{125} will be referred to throughout this paper as the shower size. \kappa scales with the curvature of the parabola which is approximately constant for all hadronic showers; thus \kappa is presently set as a default to 0.303, while the other two parameters are allowed to vary event by event. The best-fit parameters are found using a 3-step maximum-likelihood technique, which compares the timing and charge of the hits to both the expected charge LDF and an expected timing LDF. Both saturation of the tanks, as well as stations which are not hit, are taken into account in the likelihood.

The Laputop reconstruction also takes the actual snow depths on top of the tanks into account. The frozen water tanks of IceTop were deployed flush with the surface of the snow at the site. However, wind-blown snow continuously drifts over the array and covers the tanks with an overburden that increases over time and varies from tank to tank. Figure 2 shows the depths of snow covering the array’s tanks in each of the three years analyzed in this work: 2010, 2011 and 2012. The site accumulates about 20 cm of snow per year on average.

Figure 2: Depths of snow covering the IceTop tanks, measured during each of the three years studied in this work. The pink dots indicate the positions of buildings at the site. The z-axis of these plots represents the snow depth, measured in meters.

The total signal S observed by the two DOMs in the IceTop tank consists of snow attenuated electromagnetic (e^{\pm}/\gamma) particles and unattenuated muons. If the attenuation of the snow is not taken into account, the reduction in total signal amplitude will make an event look less energetic, or “smaller”, than it really is. To take snow attenuation into account, a simple exponential reduction S{}_{red} is applied to the expected S:

S_{red}=S\cdot e^{\frac{-d}{\lambda\cos\theta}}, (1)

which only depends on the slant depth of snow overburden for the tank (d/\cos\theta), and an effective attenuation length \lambda, which takes into account absorption/generation in the shower and particle-type specific absorption behavior. The same attenuation length \lambda is applied throughout the array, in the same way for all showers of any size. However, because of the increasing snow load from year-to-year, the optimal \lambda also changes from year-to-year. Thus each of the three years of data was optimized separately to find the \lambda which best creates agreement in the S_{125} spectrum across different regions in the array (deeply-buried, and sparsely-buried). These best-fit values of \lambda are: 2.1 meters for 2010/11, 2.25 meters for 2011/12, and 2.25 meters for 2012/13. Furthermore, the snow depth, d, has two sources of uncertainty De Ridder (2019). First, sastrugi333Sastrugi are irregular waves formed on the surface of the snow by wind erosion. cause variations in the snow depth across a single tank. Sastrugi heights measured in-situ were observed to follow a Gaussian distribution with a standard deviation of 4 cm. Second, the depth at each IceTop tank is measured twice per year, in February and in November, with an occasional third measurement in January; for all other times, a linear interpolation between these measurements is used to estimate the daily snow depth. Finer measurements of accumulation at the South Pole (made monthly by the Antarctic Meteorological Research Center at a site near IceTop) exhibit variations around a smooth interpolation with a \sigma of 27% of the difference between the November and February measurements (\Delta H_{snow}, in centimeters). The two sources of uncertainty are combined and the snow depth d in Eqn. 1 (the signal reduction due to the snow) is smeared by a Gaussian with a sigma of \sqrt{(4\,\text{cm})^{2}+(0.27\cdot\Delta H_{snow})^{2}} De Ridder (2019).

I.5 InIce Reconstruction

The IceCube detector observes the Cherenkov light pattern from the high-energy muon bundles that propagate through the Antarctic ice as well as their accompanying energy losses (which create Cherenkov-light-emitting cascades).

These measured observables provide a handle on the primary composition: iron-induced showers are more muon-rich than proton-induced showers of the same energy, due to the superposition model of nucleons (detailed in e.g. Kampert and Unger (2012)). Therefore, iron-induced showers are more muon-rich than proton-induced showers of the same energy, therefore iron-induced showers will result in a greater overall deposit of Cherenkov light in IceCube than proton-induced showers of the same primary energy. Furthermore, proton-induced showers are more likely than iron-induced showers to have extremely high-energy muons in the bundle, which can create larger local energy depositions from Bremsstrahlung. Consequently, proton-induced showers are expected to create fewer but higher-energy stochastic losses in the detector than iron-induced showers (of the same primary energy). On the other hand, since iron-induced showers have more total muons than those induced by protons (of the same primary energy), an iron-induced shower is likely to undergo more lower-energy stochastic losses in the detector than a proton-induced shower (of the same primary energy). Therefore, counting these stochastic fluctuations (henceforth “stochastics”) gives us additional composition-sensitive information.

In order to measure these composition-sensitive parameters for coincident events, the track position and direction as reconstructed by IceTop (with Laputop) is used, and an energy loss reconstruction algorithm (Millipede) uses the timing and charge information of the observed Cherenkov light in the ice to create a profile of energy losses along the track as a function of slant depth, with 20-meter segmentation (as discussed in detail in Aartsen et al. (2014a)).

Figure 3 shows an example energy loss profile from a cosmic ray event. Note that sections of the energy loss profile corresponding to IceCube’s dust layer444The “dust layer” is a distinctive thick layer of dust deposited several millennia ago and currently located \sim1950-2050 m beneath the surface of the ice Aartsen et al. (2013b). (where there are consequently few photons) and near the boundaries of IceCube’s volume (where reported energy losses can be irregular) are removed (as shown in gray in the figure). The gaps seen in the energy loss profile are not segments on the muon bundle track where no energy was lost, but rather are pieces of the track that are not well-sampled by the detector.

Figure 3: Example of the energy loss reconstruction of a large event, where the solid red line demonstrates the average energy loss fit, the dashed red line represents the standard stochastics selection, and the dotted red line indicates the strong stochastics selection, as noted in the legend. Note that when successive bins exceed the selection criteria they are counted individually. (For example, around slant depth of 2350 m two bins exceed the standard selection. These are therefore counted as two high energy stochastics.) The gray band is the approximate location of the dust layer for the slant depth of this particular event. Feusels (2014); De Ridder (2019)

The energy loss profile is then fit to extract two composition-sensitive parameters: a) the average energy loss behavior, which is indicated as the red line in Figure 3, and b) the size and quantity of deviations from that average behavior (the stochastics). The energy loss observable (dE_{\mu}/dX) is defined as the value of the fit to the energy loss profile at a fixed slant depth of X=1500 m, which corresponds roughly to the top of the IceCube detector (marked on the left side of Figure 3).

Two methods of selecting a number of high-energy stochastics from a energy loss profile are used in this work: a standard selection (marked as the red dashed line in Figure 3 and a strong selection requiring higher stochastic energy loss (marked as the red dotted line in Figure 3 Feusels (2014); De Ridder (2019). The selection criterion for the stochastics is given by:

\frac{\mathrm{d}E_{\mu}}{\mathrm{d}X}(X)>a\cdot\bigg{(}\frac{\mathrm{d}E_{\mu}% }{\mathrm{d}X}(X)\bigg{)}_{\mathrm{reco}}^{b},

where a=5 and b=0.8 for the standard selection, and a=7 and b=0.9 for the strong selection (with appropriate dimensions for a and b).

Figure 4: Composition sensitivity in MC simulations (Sibyll2.1) of three InIce variables: energy loss \mathrm{d}E/\mathrm{d}X (left), the number of high energy stochastics standard selection (middle), and the number of high energy stochastics strong selection (right) as discussed in Section I.5. Error bars represent RMS spread of the distribution. As the standard stochastics count begins to lose sensitivity at 100 PeV, the strong stochastics count begins to be sensitive Feusels (2014). It is clear that the energy loss is the primary composition-sensitive parameter.

Figure 4 demonstrates the composition-sensitivity of the InIce observables reconstructed by Millipede. The energy loss parameter is directly comparable to the number of high-energy muons in the air shower and is therefore the primary composition-sensitive observable, as shown in Figure 4, left. The stochastics provide additional composition-sensitivity, as shown in Figure 4, center and right: Iron bundles have a larger number of stochastic losses since they have more muons, but the energy losses from proton bundles can be more extreme since the same total energy is transferred to fewer muons. The number of stochastics with the standard selection can separate masses best at low muon multiplicities, below 100 PeV. From about 30 PeV and above, for bundles containing at least 100 muons, the stronger selection performs better. In both cases there are more high energy cascades selected for iron bundles than for proton.

The muon multiplicity in air showers detected by IceCube shows a seasonal variation, which is due to the semi-annual alternation between the polar day and night, and the accompanying temperature changes in the atmosphere. The measured variation of \log_{10}(\mathrm{d}E_{\mu}/\mathrm{d}X) is found to be 10-15% of the proton-iron difference. Since simulations are only performed with one atmosphere, the July,1997 atmosphere, all other months of data need to be corrected with respect to July. This correction correlates the changing temperature profile of the entire atmosphere, weighted with the muon production depth profile, with the measured variation of \log_{10}(\mathrm{d}E_{\mu}/\mathrm{d}X). This energy dependent correlation factor is used as a correction. A small but symmetric variation of \pm 3\% in the proton iron-space remains, which smears the data slightly. For more details, see De Ridder et al. (2013).

I.6 Quality Cuts

Quality cuts are used to ensure a sample of events which will be well reconstructed. The cuts for the IceTop-alone analysis are described in Aartsen et al. (2013a) while the cuts for the coincident analysis are detailed here.

IceTop Selections: Many of the IceTop selections revolve around the success of the Laputop reconstruction algorithm, which does an excellent job of reconstructing contained events (those with a shower core inside the area of the array) but its performance suffers for uncontained events.

  • The number of stations after cleaning is required to be \geq 5.

  • The largest snow-corrected charge measured in any tank is required to be at least 6 VEM.

  • The station with the highest deposited charge is not allowed to be at the edge of the detector.

  • The neighboring tank in the same station as the tank with the largest signal must have at least 4 VEM.

  • The fraction of hit stations within a circle centered on the center of gravity of the shower with outer radius at the furthest hit station must be greater than 0.2.

  • The reconstruction algorithm Laputop is required to converge.

  • The LDF slope parameter \beta is required to be between 1.4 and 9.5.

  • The core location of the air shower must be reconstructed within a scaled factor of 0.96 of the area of the array.

InIce Selections: The InIce reconstruction begins with a fixed track position and direction from Laputop. Therefore the InIce quality selections focus on ensuring an accurate reconstruction of the energy loss from Millipede.

  • The track position and direction calculated by Laputop is required to pass within the In-Ice instrumented volume.

  • A minimum of 8 InIce DOMs are required to be hit.

  • The Millipede energy loss reconstruction must succeed with log{}_{10}(rlogl)<2.0 and the total charge predicted (QTOT) must be at least 90% of that measured ( log{}_{10}(\frac{QTOT_{predicted}}{QTOT_{measured}})>-0.03 ) .

  • At least 3 reconstructed cascades remain after all previous selections and after removal of cascades in the dust layer and at the edge of the detector (as discussed above).

A note about the zenith angle: As in Aartsen et al. (2013a), the IceTop-alone analysis presented here is divided into four bins of zenith angle, the steepest of which is limited to cos(\theta)\geq 0.80. No explicit cut in zenith angle is applied for the coincident analysis, since the solid angle acceptance of the two detectors together limits the zenith angle range to cos(\theta)\sim 0.85.

I.7 Performance

Figure 5: Core position resolution (upper) and angular resolution (lower) of the reconstructed air shower (both defined as containing 68% of the events) as a function of the primary energy, after quality cuts.
Figure 5: Core position resolution (upper) and angular resolution (lower) of the reconstructed air shower (both defined as containing 68% of the events) as a function of the primary energy, after quality cuts.
Figure 6: Reconstructed energy loss as a function of S_{125} for proton (red) and iron (blue) simulations, with the standard deviation indicated as error bars.

The reconstruction procedure and quality cuts described above yield a set of events with a core position resolution of 6-20 meters, and a track direction resolution of 0.3-1.0 degrees. These values depend on energy: Figure 6 shows the position and angular resolutions as a function of energy, for the coincident sample and the IceTop-alone sample (divided into high-zenith and low-zenith angles). Note that the detector reaches 100% efficiency for all particle types at \sim3 PeV (Feusels (2014)). As the size of the air showers approach and exceed the size of the IceTop array and the tanks become saturated, the reconstruction becomes less precise.

Figure 6 shows the reconstructed energy loss from IceCube compared to the reconstructed S_{125}  parameter from IceTop: in this parameter space there is clearly a strong separation between proton and iron primaries.

II The IceTop-alone analysis

The IceTop-alone analysis is sensitive to the energy spectrum of cosmic rays up to the EeV energy range. The reconstructed shower size parameter (S_{125}) in particular is highly correlated to the energy of the primary air shower. Figure 7 shows the relationship between \log_{10}(S_{125}) to \log_{10}(E/\mathrm{GeV}). In Aartsen et al. (2013a), a function relating \log_{10}(S_{125}) to \log_{10}(E/\mathrm{GeV}) was derived using Monte Carlo simulations divided into four ranges of zenith angles. For a given zenith range, the distributions of true energy for each slice of 0.05 in \log_{10}(S_{125}) were weighted using the H4a model from Gaisser (2012) as a composition assumption, and fitted with a Gaussian. (Since neither silicon nor magnesium were simulated, but both are included in the H4a model, simulated oxygen was weighted by the sum of CNO and MgSi model components.) The fitted mean of the Gaussian was then used as the energy estimate for that slice in \log_{10}(S_{125}). The functional form of the conversion is:

\mathrm{log}_{10}(E/\mathrm{GeV})=p_{1}\mathrm{log}_{10}(S_{125/\mathrm{VEM}})% +p_{0}. (2)

Using updated simulation and reconstruction algorithms, the mapping of \log_{10}(S_{125}) to \log_{10}(E/\mathrm{GeV}) from Aartsen et al. (2013a) has been re-optimized and applied to the 3-year data set. The updated fit parameters for Eqn. 2 are in the following table:

 Zenith range  p_{0}  p_{1}
 0.95<\cos(\theta)\leq 1.0  6.011  0.933
 0.90<\cos(\theta)\leq 0.95  6.055  0.924
 0.85<\cos(\theta)\leq 0.90  6.110  0.915
 0.80<\cos(\theta)\leq 0.85  6.177  0.907
Table 1: Fit parameters for converting IceTop shower size S_{125}to energy in Eqn. 2, using the “H4a” composition assumption. Errors are on the order of 0.006 for p_{0}, and 0.0035 for p_{1}.
Figure 7: Relationship between S_{125} and primary energy: updated version of Figure 4 from Aartsen et al. (2013a), for primary protons at high (\cos(\theta)\geq 0.95) zenith angles.

The energy bias and resolution of this technique are shown in Figure 8. The reduced precision beyond 8.0 (100 PeV) is related to the reduced angular and position resolution shown in Figure 6, which creates an extra smearing effect in S_{125}.

Figure 8: Energy bias (upper) and resolution (lower) of the conversion functions in Table 1, for simulated events mixed according to the H4a model in the four zenith bins. (These are updated versions of Figures 7 and 8 from Aartsen et al. (2013a).)
Figure 9: All-particle energy spectrum from the IceTop-alone analysis from each of the three years, and the three years together.

Figure 9 shows the result of the IceTop-alone 3-year analysis for each year individually and combined, multiplied by a factor of E^{3} to highlight the details. The gray band represents the total systematic uncertaintity of the IceTop detector, as described in Aartsen et al. (2013a). These results are consistent between the three years and are also consistent with those previously published from one year of data Aartsen et al. (2013a).

III The Coincident analysis

When the surface observables from IceTop are combined with the additional observables from the InIce detector, the high-energy muon component of the shower is measured in coincidence with the electromagnetic component of the EAS. Using this coincident configuration, a mass-independent primary energy spectrum and individual elemental spectra are measured. This technique was developed for the measurement of the cosmic ray composition of one month of IT-40/IC-40 data in the energy range between 1 PeV to 30 PeV Andeen (2011); Abbasi et al. (2013b) using two input variables. Building on this experience, the technique was extended to five input variables over a wider primary energy range, optimized over a larger scan of different network types, and trained on more Monte Carlo simulated events. This updated technique was applied to a single year of data from the nearly complete IT-73/IC-79 detector in Feusels et al. (2013); Feusels (2014). Here, the one year analysis is improved and further expanded to include three years of data.

III.1 Neural Network Mapping Technique

This analysis includes five variables which depend on primary energy and primary mass in a non-linear fashion: the shower size in IceTop (S_{125}), the zenith angle (\cos(\theta)), the muon energy loss in the ice (\mathrm{d}E/\mathrm{d}X), and the number of high-energy stochastics under two selections (standard and strong). There is no theoretical analytical expression that relates our input variables to primary mass and primary energy; thus, an artificial neural network (NN)555In particular, a feed-forward multilayer-perceptron (MLP) neural network is used from the TMVA Hoecker et al. (2007) machine learning package. is trained on simulation to determine the relationships between the five inputs and the two outputs. The network is strongly dependent on the two primary parameters, S_{125} and \mathrm{d}E/\mathrm{d}X, but the three other parameters do contribute to the energy and mass reconstruction.

The final high-quality sample of simulated Monte Carlo data is split into three parts. Half of the sample is used to generate the neural network (the network sample). The other half (the verification sample) is used for comparisons of data and simulation in the final analysis steps. The network sample is again split in two: 74357 events are used to train the network (the training sample), the remaining 67399 events (the test sample) serve to test the network and to select the network architecture and optimal activation function based on the network performance. Networks were trained on unweighted events; however, every Monte Carlo sample mentioned above is chosen in such a way that it contains an equal mixture of each of the four primary types (p, He, O and Fe) and covers the full energy range.

During the first 5000 of 10000 minimizer iterations (also called cycles or epochs), only a random selection of 60% of the training data is utilized. After the training converged on this random selection, the training continues on the full training set.

III.2 Optimizing the Neural Network

Many different neural network architectures were evaluated for performance before analyzing any data, as discussed in Feusels (2014). In addition to networks with 5 inputs as described above, alternative networks with the 2 primary inputs (\log_{10}(S_{125}) and \log_{10}(\mathrm{d}E/\mathrm{d}X) only), 3 inputs (adding \cos(\theta)) and 4 inputs (adding the standard selection of high energy stochastics only) were tested. Three groups of network structures were explored: with one, two, and three hidden layers, and the number of neurons was varied within the hidden layers. Two activation functions (a sigmoid, and a tanh) were explored. In total, 207 networks for each of the two activation functions and for each number of inputs (1656 networks in total) were trained on the simulations.

The performance of each network was assessed according to how well it reconstructed primary energy and primary mass. The assessment process was optimized to find the network with the smallest and most consistent RMS spread and bias over all energies, and which had mass groups that were best-separated and most distinctive (i.e. “peaky”). The final optimized network has 5 inputs, 7 neurons in a first hidden layer, 4 neurons in a second hidden layer, and 2 outputs, with a tanh activation function connecting the neurons and a linear mapping from the last layer to the output neurons. A schematic of this network is shown in Figure 10.

Figure 10: The neural network architecture of the best performing neural network. This network maps five input variables onto two output variables using two hidden layers with respectively seven and four neurons using a tanh activation function. It is therefore called a 5-7-4-2 network.

It is important to note that this neural network has two target outputs which are very different in nature: the first output is a continuous energy distribution, the second target output is instead is composed of four discrete numbers corresponding to four elemental masses simulated. Therefore, the neural network energy output (E{}_{0,reco}) is also a continuous distribution which is expected to reproduce the true primary energy (within some bias and resolution) for each event, as discussed below in Sec. III.3. On the other hand, the neural network mass output results in smeared distributions around the four discrete mass numbers, which require further analysis in order to decompose the primary mass. The mass is therefore not reconstructed on an event-by-event basis but is determined statistically for the entire data set, as discussed below in Sec. III.4.

III.3 Neural Network Primary Energy Reconstruction

Figure 11: Energy reconstruction bias (upper) and resolution (lower) as a function of the reconstructed energy for the different primary types and for an equal mixture of each type.

The energy dependence of the primary energy bias and resolution as reconstructed by the NN are shown in Figure 11. The energy resolution (Figure 11, lower) ranges from 9% (for iron showers at around 30 PeV) and 18% with the worst resolutions below the energy threshold of \sim3 PeV and at the highest energies due to the worsening core position and angular resolution (as discussed in Section I.7). Heavier primaries can be reconstructed more precisely because of their lower intrinsic shower fluctuations. As mentioned in Section II, the overall decrease in precision beyond \sim100 PeV is partially caused by the decrease in precision in angular and position resolution shown in Figure 6, which creates an extra smearing effect in S_{125}.

In this analysis, events are divided into energy bins of width 0.1 in \log_{10}(E/\mathrm{GeV}) , which is larger than both the energy bias and the energy resolution as shown in Figure 11. However, due to the decrease in accuracy, precision, and available statistics at high energies (\log_{10}(E/\mathrm{GeV}) > 8.0), bins of width 0.2 are used in this region. Above 1 EeV the energy bias dependence on the primary type becomes too large and limits the energy range over which this analysis is optimal.

Figure 12 shows the all-particle energy spectrum results for the coincident analysis for the three years individually and combined, multiplied by a factor of E^{3} to highlight the details: the results are consistent between the years. The gray band represents the combined systematic uncertainties of the IceTop and InIce detectors for the coincident analysis, as discussed in Section IV.2. These results are included in Table 3 in Appendix .4.

Figure 12: All-particle energy spectrum from the coincident analysis from each of the three years analyzed individually compared to the combined result. The gray band represents the total detector uncertainty from both the IceTop and InIce arrays, as discussed in Section IV.2.

III.4 Composition Reconstruction using kernel density estimation to fit neural network templates

Figure 13 shows histograms for each simulated element (proton, helium, oxygen and iron) in the natural logarithm of the neural network mass output for one slice in reconstructed energy. (The four simulated types (proton, helium, oxygen and iron) are equidistant in \langle\ln A\rangle, but not in A. Thus, the histograms are expected to be more distinct in logarithmic space.) In every slice in energy, the histogram for each primary element is converted into a template probability density function (p.d.f.) using an adaptive kernel density estimation (KDE) method Cranmer (2001). The template p.d.f.’s are shown as the solid lines in Figure 13. The template p.d.f.’s for all energy slices used in this analysis are given in Appendix .5 in Figure 25. The four primary types exhibit four very distinctive shapes in each slice in energy over the whole energy range. At \log_{10}(E/\mathrm{GeV}) = 8.0, the template p.d.f.’s begin to exhibit greater overlap due to the limited statistics in the Monte Carlo sets, reducing the composition sensitivity of the analysis. Beyond \log_{10}(E/\mathrm{GeV}) = 9.0 (1 EeV), the analysis becomes unreliable due to the overlap and reduction in data statistics.

Figure 13: Example distribution of the natural logarithm of the neural network mass output in the energy range between 7.4 and 7.5 in log{}_{10}(E{}_{0,reco}/GeV). The y-axis represents the number of simulated events for proton (red), helium (orange), oxygen (green) and iron (blue). The solid lines represent the probability density function (p.d.f.) found by the adaptive KDE method.

Using the Roofit package Verkerke and Kirkby (2003), the set of four template p.d.f.’s were then weighted to find the fractions (which were constrained to add to unity) which best fit the NN mass output for the experimental data in the same slice in reconstructed energy666This serves a similar purpose to the chi-squared minimization approach described in Andeen (2011); Abbasi et al. (2013b); however, the new unbinned extended likelihood technique improves on the previous method by correctly taking into account the Poisson fluctuations of the bin contents in both the data and the templates, which is particularly relevant in bins with few events.. The result of this method applied to the experimental data for the same energy bin as in Figure 13 is shown in Figure 14 (upper). Additionally, the correlation between the fitted weights is shown in the form of uncertainty elipses in Figure 14 (lower). (The fit plots for all energy bins are given in Figure 26 in Appendix .5 and the corresponding correlation coefficients are shown in Table 6 in Appendix .5.) The resulting fractions of neighboring elements (i.e. protons and helium) are anti-correlated, while those from distant elements (i.e. protons and iron) are virtually uncorrelated: this means a proton primary is more likely to be confused for a helium primary than for an iron primary, which is expected.

Figure 14: upper: Example data distribution fit of the natural logarithm of the neural network mass output in the energy range between 7.4 and 7.5 in log{}_{10}(E{}_{0,reco}/GeV) and lower: corresponding contour plot with the best fit values. Data distributions, contour plots, and complete correlation coefficient matrices for all energy bins are included in Appendix .5.

The KDE template-fitting procedure yields a measurement of the fractions of each of the four nuclear mass groups (represented by H, He, O, and Fe), for each bin in energy. Each of these four individual fractions is then translated into an individual spectrum for the corresponding elemental group, as shown in Figure 15 (colors) compared to the all-particle spectrum (black). Recent model predictions are also included in Figure 15, which will be discussed in Section V. The gray band represents the total coincident detector uncertainty from both the IceTop and InIce arrays, which will be discussed in Section IV.2. These results are included in Tables 4-5 in Appendix .4.

Intermediate elements, not part of the four groups listed above, are expected to produce neural network outputs in between the adjacent groups, so will partially contribute to the flux of the groups that bracket it. In order to test this, a small sample of silicon was passed through the NN + KDE chain and treated as “data”. The natural log of the mass of silicon is approximately midway between that of oxygen and that of iron; therefore, as expected (due to the regression-style neural network mass output), the silicon is reconstructed as a nearly 50/50 mixture of oxygen and iron across all energies.

Figure 15: Individual spectra for the four mass groups (protons in red, helium in yellow, oxygen in green, and iron in blue) including total detector systematic compared with various cosmic ray models (H3a and H4a Gaisser (2012)) and phenomenological experimental fits (GSTGaisser et al. (2013) and GSF Dembinski et al. (2018)).

Figure 16 shows the mean log mass, which is derived from the individual fractions shown in Figure 15. Again, the gray band represents the total coincident detector uncertainty from both the IceTop and InIce arrays, which will be discussed in Section IV.2. Each of the three years of data are again shown both separately and combined and agree very well within the statistical and systematic uncertainties.

Figure 16: Mean log mass (\langle\ln A\rangle) for the three individual years, and the three years combined. The gray band represents the combined systematic uncertainties of the IceTop and InIce detectors for the coincident analysis, as discussed in Section IV.2.

IV Systematic Uncertainties

The uncertainties in the coincident analysis reported here can be grouped into three categories: analysis method, detector effects, and the hadronic interaction model.

IV.1 Analysis Method

As described above, the shape of the Monte Carlo templates is derived from an adaptive KDE method, which determines the optimal width of the Gaussian kernel function. To check the robustness of the composition fitting results, the optimal kernel width is artificially modified by a factor of \pm 90\%, resulting in either very jagged or very smooth templates. These artificially modified templates are then used in place of the optimal templates for the remainder of the analysis in order to measure the effect of the vastly different templates on the results. The variation in the final results due to the modified template shapes is so small that it is not visible in a figure; however, the values for uncertainty are included in the tables in Appendix .5.

IV.2 Detector Uncertainties

Three main detector effects contribute to the uncertainty in the composition results: the snow correction, the absolute energy scale of IceTop, and the light yield in the ice.

Figure 17: Total detector uncertainty (gray) on all-particle energy spectrum from the combination of light yield (magenta), snow correction (cyan) and energy scale uncertainty (orange).
Figure 18: Total detector uncertainty (gray) on \langle\ln A\rangle from the combination of light yield (magenta), snow correction (cyan) and energy scale uncertainty (orange).
Figure 19: Total detector uncertainty (gray) on elementary energy spectra from the combination of light yield (magenta), snow correction (cyan) and energy scale uncertainty (orange).

IV.2.1 Snow correction

Although the method for snow correction described in Section I.4 works well on average over the entire energy and zenith angle region, it is not perfect. This is predominantly due to its inability to distinguish between the electromagnetic and muonic component of the air shower. A systematic uncertainty of \pm0.2 m was assigned to \lambda, which covers the variations due to the unknown composition, energy and zenith angle dependence, as discussed in the Appendix of Aartsen et al. (2013a).

IV.2.2 Absolute energy scale

To obtain a reliable energy scale, the simulated response of each IceTop tank is calibrated using the signal from atmospheric muons Abbasi et al. (2013a). This procedure has proven to be very reliable, but due to the unknown composition, atmospheric conditions, etc, a final 3% uncertainty on the charge calibration, and thus on the absolute energy scale, needs to be taken into accountVan Overloop (2011). This translates directly to a 3% shift in S{}_{125}.

IV.2.3 Light Yield

The uncertainty on the photon detection efficiency by the DOM (referred to as the DOM efficiency) is \pm3%, which implies a \pm3% possible variation of the total light yield observed. The uncertainty caused by the photon propagation in the ice includes scattering and absorption coefficients in the bulk ice and an effective scattering length of the hole ice, the ice in the drill hole around the DOM Rongen et al. (2016). For the bulk ice scattering and absorption coefficient, 3 points are taken on the 1 \sigma error ellipse around the nominal values, shown in Aartsen et al. (2013c): a +10% scattering coefficient, a +10% absorption coefficient and a -7.1% scattering and absorption coefficient. Alternate effective scattering lengths of 30 cm and 100 cm were used for the hole ice model. For all these individual systematic uncertainties, simulations were produced and their effect was studied. It was found that the main combined effect is to influence the light yield in the DOM, and that this effect is rather independent of the initial light yield and zenith angle of the muon bundle. This means that all those systematic uncertainties can be combined and modeled as a shift on the observed light yield. Furthermore, the systematic uncertainties are (nearly) uncorrelated; thus the various errors and shift on light yield are added in quadrature, which gives a total light yield uncertainty of +9.6% and -12.5%. The individual contribution to the observed light yield shifts, as well as the total light yield uncertainty are given in Table 2.

Table 2: Sytematic light yield shift
Effect Light yield shift
+10% scattering +3.6%
+10% absorption -11.8%
-7.1% scattering and absorption +7%
30 cm hole ice scattering +4.5%
100 cm hole ice scattering -2.9%
DOM efficiency \pm 3\%
Total Light Yield Effect +9.6%,-12.5%

IV.2.4 Total Detector Uncertainties

The three detector uncertainties discussed above (snow correction, absolute energy scale, and light yield) are added together in quadrature into a total detector uncertainty. Figures 1718, and 19 show individual and combined contributions to the uncertainty in the all-particle energy spectra, mean log mass, and the individual elemental spectra. In all figures, the gray bands are the total combined detector uncertainty which match the gray bands shown in Figures 12, 16, and 15, respectively, and are included in the tables in Appendix .4.

IV.3 Hadronic Interaction Model

The influence of the hadronic interaction model on the measurement of the cosmic ray composition arises mainly from variations between the models in the predicted number of high-energy muons. Sibyll 2.1 was used as the hadronic interaction model for the baseline simulation datasets, while the three post-LHC models–QGSJET II-04Ostapchenko (2011), Sibyll 2.3Riehn et al. (2016) and EPOS-LHCPierog et al. (2015)–are used as alternate models. The variation between the models in the high-energy (>300 GeV) muon number is smaller than for the surface (GeV) muons, but is still of order 15% at maximum (between QGSJet II-04 and EPOS-LHC).

It is important to note that the main Sibyll 2.1 dataset included full samples across all energies for each particle type (proton, helium, oxygen, and iron); however, due to the computational time involved in generating a full simulated data set (\sim60000 CPU years), the alternate models include 1/10th the number of events, in proton and iron only. As a result of the limited size of the alternate simulated data sets, a full analysis is not repeated using those datasets for the simulated templates. Rather, the difference in the number of high-energy muons and in S_{125} in the alternative sample with respect to Sibyll 2.1 is calculated for each particle type and applied to the experimental data, which is then passed through the full analysis chain. (Neither the stochastics nor the zenith angle are taken into account in this estimate due to their low impact on the results.) The results of this process are then weighted using the reconstructed elemental fractions for each particle type from the baseline Sibyll 2.1 analysis in order to obtain an estimate of what the final result would have been had the alternative interaction model been used as the nominal Monte Carlo simulation instead of Sibyll 2.1. It is probable that the alternate results reported here differ slightly from the actual fractions that would be reconstructed with a full alternative simulated dataset.

The estimated uncertainty due to each alternative model is shown in Figures 20, 21, and 22, which depict the effect of the choice of hadronic interaction model in the all-particle energy spectra, mean log mass, and the individual elemental spectra.

Figure 20: Hadronic interaction model uncertainty range on all-particle energy spectrum based on EposLHC (blue), Sibyll2.3 (red) and QGSJetII-04 (green).
Figure 21: Hadronic interaction model uncertainty range on \langle\ln A\rangle  based on EposLHC (blue), Sibyll2.3 (red) and QGSJetII-04 (green).
Figure 22: Hadronic interaction model uncertainty range on elementary energy spectra based on EposLHC (blue), Sibyll2.3 (red) and QGSJetII-04 (green).

V Results, Discussion and Outlook

V.1 Energy Spectra

Figure 23 compares the energy spectra resulting from the IceTop-alone analysis and the coincident analysis as described herein. The analyses are consistent with each other within the statistical and systematic uncertainties (only the smaller IceTop-alone systematic uncertainties are shown, for clarity). This good agreement indicates that the dependence of the IceTop-alone analysis on composition model has been effectively mitigated through the analysis technique, since the IceTop/InIce coincident analysis doesn’t require prior knowledge of the composition.

Figure 23: A comparison of the combined three-year spectra from the two analyses in this paper: the IceTop-alone analysis, and the coincident analysis. The gray band represents the total systematic uncertaintity of the IceTop detector, as described in Aartsen et al. (2013a)

V.2 Mass composition

Figure 24: Comparison of the all-particle and composition spectra of the four elemental groups H, He, O and Fe from this analysis using Sibyll 2.1 (black) with other experiments. The data set for the all-particle spectra are taken from Tibet Amenomori et al. (2008), Tunka Prosin, V.V. et al. (2016), Yakutsk Ivanov et al. (2009), Tale Abbasi et al. (2018), Hires Abbasi et al. (2008b) and Telescope Array Abbasi et al. (2016). The Kascade Antoni et al. (2005) results are using a 5 component fit of H, He, CNO, MgSi and Fe groups using Sibyll 2.1. Therefore only the H and He spectra are compared directly as the other groups are strongly correlated. Kascade-Grande Arteaga-Velázquez et al. (2018) results are using a 3 component fit of H, HeCNO and Heavy groups using Sibyll 2.3. Therefore only the H and Heavy spectra are compared, as the HeCNO group cannot be deconvoluted into a He and CNO part. The Tunka Epimakhov et al. (2013) results are using a 4 component fit of H, He, N and Fe groups using QGSJET II-04. The Pierre Auger Observatory Bellido et al. (2018) results are calculated by using their published elementary group fraction for H, He, N and Fe using Sibyll 2.3 convoluted with their most recent energy spectrum. Note that differences in how different experiments handle intermediate elements (not one of the four groups used here) may lead to some small systematic differences in flux measurements between different experiments.

The elemental spectra results, shown in Figure 15, agree well with the recent H3a and H4a phenomenological models of the transition region between galactic and extragalactic cosmic rays (Gaisser (2012)), in which heavier elements retain a harder spectral index to higher energies. Within the statistical and systematical uncertainties the elemental spectra are also compatible with the phenomenological data fits, Gaisser-Stanev-Tilav (GST) Gaisser et al. (2013) fit and Global Spline Fit (GSF) Dembinski et al. (2018). These elemental spectra correlate with an increase in the mean-log-mass as a function of energy until about 100–200 PeV, as shown in Figure 16 (which is also derived from the individual fractions). Beyond this energy, given the statistical and systematic uncertainties the data are consistent with a composition that is either unchanging or decreasing. In Figure 24, our results are compared with those from other recent experiments: the results reported here indicate a higher flux in the iron group at high energies. As shown in Section IV.3, the absolute scale of the composition is strongly dependent on which hadronic interaction model is used for the simulations. In fact, the uncertainty due to the choice of hadronic interaction model is the biggest limitation on our analysis.

V.3 Outlook

The all-particle energy spectrum results presented here are consistent with each other and with previously published IceCube results Aartsen et al. (2013a). The limiting systematic effect is the uncertainty of the snow coverage over the tanks. The composition analysis results presented here are significantly improved from previously published results, which included only one month of data taken with a partly completed array Andeen (2011); Abbasi et al. (2013b); however, the present results are still limited by the amount of data on hand, the systematic uncertainty due to detector effects (particularly the light yield in the ice), and the dependence on the choice of hadronic interaction model used for the simulations. For future analyses, we plan to include more years of experimental data, to simulate more intermediate elements, to investigate new composition-sensitive parameters currently under development, and to incorporate results from new internal studies to reduce the detector systematic uncertainties. These updates will improve the precision of both analyses, and enable the extension of the analyses to higher and lower energies. Furthermore, the analyses presented here are well-suited to capitalize on future extensions to the IceCube Neutrino Observatory Aartsen et al. (2014b).

Acknowledgements.
The IceCube collaboration acknowledges the significant contributions to this manuscript from Marquette University, the University of Alaska Anchorage, and the University of Gent. This paper also profited enormously from the contributions of our departed colleague Stefan Westerhoff (1967-2018). USA – U.S. National Science Foundation-Office of Polar Programs, U.S. National Science Foundation-Physics Division, Wisconsin Alumni Research Foundation, Center for High Throughput Computing (CHTC) at the University of Wisconsin-Madison, Open Science Grid (OSG), Extreme Science and Engineering Discovery Environment (XSEDE), U.S. Department of Energy-National Energy Research Scientific Computing Center, Particle astrophysics research computing center at the University of Maryland, Institute for Cyber-Enabled Research at Michigan State University, and Astroparticle physics computational facility at Marquette University; Belgium – Funds for Scientific Research (FRS-FNRS and FWO), FWO Odysseus and Big Science programmes, and Belgian Federal Science Policy Office (Belspo); Germany – Bundesministerium für Bildung und Forschung (BMBF), Deutsche Forschungsgemeinschaft (DFG), Helmholtz Alliance for Astroparticle Physics (HAP), Initiative and Networking Fund of the Helmholtz Association, Deutsches Elektronen Synchrotron (DESY), and High Performance Computing cluster of the RWTH Aachen; Sweden – Swedish Research Council, Swedish Polar Research Secretariat, Swedish National Infrastructure for Computing (SNIC), and Knut and Alice Wallenberg Foundation; Australia – Australian Research Council; Canada – Natural Sciences and Engineering Research Council of Canada, Calcul Québec, Compute Ontario, Canada Foundation for Innovation, WestGrid, and Compute Canada; Denmark – Villum Fonden, Danish National Research Foundation (DNRF), Carlsberg Foundation; New Zealand – Marsden Fund; Japan – Japan Society for Promotion of Science (JSPS) and Institute for Global Prominent Research (IGPR) of Chiba University; Korea – National Research Foundation of Korea (NRF); Switzerland – Swiss National Science Foundation (SNSF).

References

  • Aartsen et al. (2017) M. Aartsen et al. (IceCube), J. Inst. 12, P03012 (2017).
  • Abbasi et al. (2008a) R. Abbasi et al. (IceCube), Nucl. Instr. and Meth. A 601, 294 (2008a).
  • Abbasi et al. (2009) R. Abbasi et al. (IceCube), Nucl. Instr. and Meth. A 618, 139 (2009).
  • Abbasi et al. (2013a) R. Abbasi et al. (IceCube), Nucl. Instr. and Meth. A 700, 188 (2013a).
  • Aartsen et al. (2013a) M. Aartsen et al. (IceCube), Physical Review D 88, 042004 (2013a).
  • Abbasi et al. (2013b) R. Abbasi et al. (IceCube), Astroparticle Physics 42, 33 (2013b).
  • Feusels et al. (2013) T. Feusels, K. Rawlins, et al.Proceedings, 33{}^{rd} International Cosmic Ray Conference ICRC2013, 861 (2013).
  • Feusels (2014) T. Feusels, Measurement of cosmic ray composition and energy spectrum between 1 PeV and 1 EeV with IceTop and IceCube, Ph.D. thesis, Gent Uni. (2014).
  • Heck et al. (1998) D. Heck et al.CORSIKA: A Monte Carlo code to simulate extensive air showers, Report FZKA 6019, Forschungszentrum Karlsruhe (1998).
  • (10) D. Heck and T. Pierog, Extensive air shower simulation with CORSIKA: A user’s guide.
  • Hillas (1997) A. M. Hillas, Nucl. Phys. Proc. Suppl. 52B, 29 (1997).
  • Battistoni et al. (2007) G. Battistoni et al., AIP Conference Proceedings 896, 31 (2007).
  • Ahn et al. (2009) E. Ahn, R. Engel, T. Gaisser, P. Lipari,  and T. Stanev, Physical Review D 80, 94003 (2009).
  • Nelson et al. (1985) W. Nelson, H. Hirayama,  and D. Rogers, report SLAC-265  (1985).
  • Hedin (1991) A. E. Hedin, Journal of Geophysical Research 96, 1159 (1991).
  • Agostinelli et al. (2003) S. Agostinelli et al., Nucl. Instr. and Meth. A 506, 250 (2003).
  • Allison et al. (2006) J. Allison et al., IEEE Transactions on Nuclear Science 53, 270 (2006).
  • Chirkin and Rhode (2004) D. Chirkin and W. Rhode, preprint  (2004), hep-ph/0407075.
  • Lundberg (2007) J. Lundberg, Nucl. Instr. and Meth. A A581 (2007).
  • Aartsen et al. (2013b) M. G. Aartsen et al. (IceCube), Nucl. Instr. and Meth. A 711, 73 (2013b).
  • Klepser (2008) S. Klepser, Reconstruction of Extensive Air Showers and Measurement of the Cosmic Ray Energy Spectrum in the Range of 1-80 PeV at the South Pole, Ph.D. thesis, Humbolt-Universität zu Berlin (2008).
  • De Ridder (2019) S. De Ridder, Sensitivity of IceCube Cosmic Ray measurements to the hadronic interaction models, Ph.D. thesis, Gent Uni. (2019).
  • Kampert and Unger (2012) K.-H. Kampert and M. Unger, Astroparticle Physics 35, 660 (2012).
  • Aartsen et al. (2014a) M. G. Aartsen et al. (IceCube), J. Inst. 9, P03009 (2014a).
  • De Ridder et al. (2013) S. De Ridder, T. Feusels, et al.Proceedings, 33{}^{rd} International Cosmic Ray Conference ICRC2013, 763 (2013).
  • Gaisser (2012) T. K. Gaisser, Astroparticle Physics 35, 801 (2012).
  • Andeen (2011) K. G. Andeen, First Measurements of Cosmic Ray Composition from 1-50 PeV Using New Techniques on Coincident Data from the IceCube Neutrino Observatory, Ph.D. thesis, University of Wisconsin-Madison (2011).
  • Hoecker et al. (2007) A. Hoecker, P. Speckmayer, J. Stelzer, J. Therhaag, E. von Toerne,  and H. Voss, PoS ACAT, 040 (2007).
  • Cranmer (2001) K. Cranmer, Computer Physics Communications 136, 198 (2001).
  • Verkerke and Kirkby (2003) W. Verkerke and D. P. Kirkby, eConf C0303241, MOLT007 (2003).
  • Gaisser et al. (2013) T. K. Gaisser, T. Stanev,  and S. Tilav, Front. Phys. (Beijing) 8, 748 (2013).
  • Dembinski et al. (2018) H. P. Dembinski, R. Engel, A. Fedynitch, T. Gaisser, F. Riehn,  and T. Stanev, PoS ICRC2017, 533 (2018).
  • Van Overloop (2011) A. Van Overloop, Proceedings, 32{}^{nd} International Cosmic Ray Conference ICRC2011, 97 (2011).
  • Rongen et al. (2016) M. Rongen et al.EPJ Web of Conferences 116, 06011 (2016).
  • Aartsen et al. (2013c) M. G. Aartsen et al. (IceCube), Nucl. Instr. and Meth. A 711, 73 (2013c).
  • Ostapchenko (2011) S. Ostapchenko, Phys. Rev. D 83, 014018 (2011).
  • Riehn et al. (2016) F. Riehn, R. Engel, A. Fedynitch, T. K. Gaisser,  and T. Stanev, PoS ICRC2015, 558 (2016).
  • Pierog et al. (2015) T. Pierog, I. Karpenko, J. M. Katzy, E. Yatsenko,  and K. Werner, Phys. Rev. C 92, 034906 (2015).
  • Amenomori et al. (2008) M. Amenomori et al.The Astrophysical Journal 678, 1165 (2008).
  • Prosin, V.V. et al. (2016) Prosin, V.V. et al.EPJ Web of Conferences 121, 03004 (2016).
  • Ivanov et al. (2009) A. A. Ivanov et al.New Journal of Physics 11, 065008 (2009).
  • Abbasi et al. (2018) R. U. Abbasi et al. (Telescope Array), The Astrophysical Journal 865, 74 (2018).
  • Abbasi et al. (2008b) R. U. Abbasi et al. (High Res.), Phys. Rev. Lett. 100, 101101 (2008b).
  • Abbasi et al. (2016) R. Abbasi et al. (Telescope Array), Astroparticle Physics 80, 131 (2016).
  • Antoni et al. (2005) T. Antoni et al. (KASCADE), Astroparticle Physics 24, 1 (2005).
  • Arteaga-Velázquez et al. (2018) C. J. Arteaga-Velázquez et al.PoS ICRC2017, 316 (2018).
  • Epimakhov et al. (2013) S. Epimakhov et al.Proceedings, 33{}^{rd} International Cosmic Ray Conference ICRC2013, 326 (2013).
  • Bellido et al. (2018) J. Bellido et al.PoS ICRC2017, 506 (2018).
  • Aartsen et al. (2014b) M. G. Aartsen et al. (IceCube),  (2014b), arXiv:1412.5106 [astro-ph.HE] .

.4 Table of Results

The values of the all particle energy spectrum from the coincident analysis, and the corresponding statistical and detector-related systematic uncertainties, as shown in Figure 12, are listed in Table 3. The values for the elemental fluxes, as shown in Figure 15, are listed in Tables 4 and 5.

Table 3: Total flux.
Energy Range Flux \pm Stat. + Det.Syst. - Det.Syst.
(log{}_{10}(E{}_{\mathrm{reco}} /GeV)) GeV{}^{-1} m{}^{-2} s{}^{-1} sr{}^{-1} GeV{}^{-1} m{}^{-2} s{}^{-1} sr{}^{-1} GeV{}^{-1} m{}^{-2} s{}^{-1} sr{}^{-1} GeV{}^{-1} m{}^{-2} s{}^{-1} sr{}^{-1}
6.55 9.0049e-14 8.4237e-17 5.8687e-15 6.6752e-15
6.65 4.6120e-14 5.3830e-17 3.1722e-15 3.7672e-15
6.75 2.3234e-14 3.4246e-17 1.6810e-15 2.0252e-15
6.85 1.1584e-14 2.1716e-17 8.7362e-16 1.0305e-15
6.95 5.6351e-15 1.3548e-17 4.4578e-16 5.2388e-16
7.05 2.6977e-15 8.3460e-18 2.1966e-16 2.5376e-16
7.15 1.3083e-15 5.1748e-18 9.4404e-17 1.1737e-16
7.25 6.4728e-16 3.2409e-18 4.7338e-17 5.9226e-17
7.35 3.2496e-16 2.0445e-18 2.3457e-17 2.6609e-17
7.45 1.6572e-16 1.3000e-18 1.1554e-17 1.3341e-17
7.55 8.6003e-17 8.3382e-19 5.9613e-18 6.8521e-18
7.65 4.3990e-17 5.3096e-19 3.7695e-18 4.3033e-18
7.75 2.1944e-17 3.3389e-19 1.3697e-18 2.1877e-18
7.85 1.1324e-17 2.1356e-19 9.8301e-19 8.6129e-19
7.95 5.5318e-18 1.3290e-19 4.8324e-19 5.7379e-19
8.10 2.0619e-18 4.8202e-20 1.4835e-19 2.0509e-19
8.30 4.6190e-19 1.8086e-20 4.2901e-20 5.4912e-20
8.50 9.9658e-20 6.6602e-21 1.0195e-20 1.4952e-20
8.70 2.2371e-20 2.5016e-21 1.2386e-21 2.0616e-21
8.90 5.2831e-21 9.6379e-22 7.2014e-22 1.0212e-21
Table 4: Proton and helium group flux.
Energy Range Flux + Stat. - + Det.Syst. -
(log{}_{10}(E{}_{\mathrm{reco}} /GeV)) GeV{}^{-1} m{}^{-2} s{}^{-1} sr{}^{-1} GeV{}^{-1} m{}^{-2} s{}^{-1} sr{}^{-1} GeV{}^{-1} m{}^{-2} s{}^{-1} sr{}^{-1}
Proton
6.55 2.7906e-14 4.4786e-15 4.7067e-15 1.2534e-14 1.2888e-14
6.65 1.3496e-14 1.6966e-15 1.7348e-15 6.3740e-15 6.2854e-15
6.75 7.8379e-15 1.0555e-15 1.0919e-15 3.8338e-15 3.7944e-15
6.85 3.4563e-15 5.4215e-16 5.6256e-16 2.0130e-15 2.0453e-15
6.95 1.3932e-15 2.0142e-16 2.0666e-16 7.8527e-16 8.1643e-16
7.05 6.6844e-16 1.4620e-16 1.5496e-16 4.1056e-16 4.3345e-16
7.15 2.6561e-16 6.0132e-17 6.3519e-17 1.7713e-16 1.9418e-16
7.25 1.0401e-16 2.2997e-17 2.3054e-17 6.4190e-17 7.2726e-17
7.35 5.2912e-17 1.5563e-17 1.6705e-17 3.6871e-17 4.3973e-17
7.45 3.1722e-17 8.0669e-18 8.5088e-18 2.2176e-17 2.1677e-17
7.55 1.2937e-17 4.0482e-18 4.1211e-18 1.0820e-17 1.0349e-17
7.65 8.2046e-18 1.6089e-18 1.0219e-18 4.8378e-18 4.4545e-18
7.75 2.5361e-18 6.4621e-19 6.9071e-19 1.1671e-18 1.8659e-18
7.85 1.1937e-18 3.7941e-19 4.3092e-19 7.9555e-19 1.2367e-18
7.95 8.3181e-19 2.9567e-19 2.5739e-19 6.2558e-19 2.8144e-19
8.10 2.4217e-19 1.1176e-19 7.0320e-20 1.4412e-19 1.1033e-19
8.30 5.2570e-20 1.7814e-20 1.5887e-20 3.6627e-20 1.8596e-20
8.50 1.3401e-20 4.3691e-21 4.0624e-21 6.3068e-21 8.5046e-21
8.70 3.5617e-21 2.8248e-21 3.3161e-21 2.9081e-21 5.4622e-21
8.90 1.6361e-23 1.6361e-23 2.3193e-21 0.0000e+00 2.1155e-21

Helium
6.55 2.2320e-14 9.8057e-15 9.7445e-15 7.5798e-15 5.1364e-15
6.65 1.5202e-14 4.0093e-15 4.1164e-15 3.7262e-15 3.6669e-15
6.75 4.7981e-15 2.1095e-15 2.1679e-15 2.1518e-15 1.9837e-15
6.85 2.8014e-15 1.1658e-15 1.1614e-15 1.3592e-15 1.1742e-15
6.95 1.6819e-15 4.0626e-16 4.1180e-16 3.0194e-16 2.1932e-16
7.05 6.0673e-16 3.1203e-16 3.0493e-16 2.2788e-16 1.5485e-16
7.15 4.1458e-16 1.4075e-16 1.3990e-16 7.9148e-17 5.9167e-17
7.25 2.2406e-16 6.4867e-17 6.7408e-17 2.8992e-17 1.0523e-17
7.35 7.8931e-17 3.6026e-17 3.5255e-17 2.5378e-17 1.9071e-17
7.45 2.1944e-17 1.6932e-17 1.6684e-17 7.3765e-18 1.1212e-17
7.55 1.4821e-17 7.0807e-18 7.1988e-18 5.1506e-18 7.1775e-18
7.65 1.1842e-19 1.1842e-19 3.1426e-18 1.1842e-19 2.0605e-18
7.75 2.3156e-18 1.3796e-18 1.3792e-18 9.4788e-19 4.4872e-19
7.85 1.4020e-18 9.2161e-19 8.5754e-19 1.0573e-18 1.0964e-19
7.95 3.1820e-19 3.1820e-19 6.3070e-19 0.0000e+00 8.3918e-19
8.10 8.1497e-21 8.1497e-21 2.7743e-19 8.1497e-21 2.1588e-19
8.30 5.0014e-25 5.0014e-25 3.2636e-20 0.0000e+00 7.4760e-20
8.50 1.1816e-24 1.1816e-24 7.6326e-21 1.1816e-24 6.9811e-23
8.70 4.0957e-21 4.0957e-21 4.0937e-21 4.0957e-21 2.1715e-21
8.90 2.7852e-21 2.7852e-21 1.1812e-21 2.7852e-21 1.6489e-22
Table 5: Oxygen and iron group fluxes.
Energy Range Flux + Stat. - + Det.Syst. -
(log{}_{10}(E{}_{\mathrm{reco}} /GeV)) GeV{}^{-1} m{}^{-2} s{}^{-1} sr{}^{-1} GeV{}^{-1} m{}^{-2} s{}^{-1} sr{}^{-1} GeV{}^{-1} m{}^{-2} s{}^{-1} sr{}^{-1}
Oxygen
6.55 1.9747e-14 9.7185e-15 9.4212e-15 8.5661e-15 4.7599e-15
6.65 6.7138e-15 5.3052e-15 5.0970e-15 5.9089e-15 2.5963e-15
6.75 5.8033e-15 2.0678e-15 1.9783e-15 2.0396e-15 1.0165e-15
6.85 2.7657e-15 1.1149e-15 1.1081e-15 1.1372e-15 6.8890e-16
6.95 1.2059e-15 3.9459e-16 3.8636e-16 3.5492e-16 1.1224e-16
7.05 8.7568e-16 3.0373e-16 3.0551e-16 2.8525e-16 1.5737e-16
7.15 2.8307e-16 1.3536e-16 1.3522e-16 1.5314e-16 5.8205e-17
7.25 9.0022e-17 7.8710e-17 7.6269e-17 9.0022e-17 6.1959e-17
7.35 9.9338e-17 3.6632e-17 3.7124e-17 7.0244e-17 3.5200e-17
7.45 6.7229e-17 1.5277e-17 1.5375e-17 3.4903e-17 1.4413e-17
7.55 2.6472e-17 6.2532e-18 6.1411e-18 2.0952e-17 1.0616e-17
7.65 2.1996e-17 3.1135e-18 2.0254e-18 1.1487e-17 3.4055e-18
7.75 9.1884e-18 1.4926e-18 1.5172e-18 4.5858e-18 3.1836e-18
7.85 3.8662e-18 8.4914e-19 8.8871e-19 2.9245e-18 2.0327e-18
7.95 1.7920e-18 6.0110e-19 5.2246e-19 1.7920e-18 7.1089e-19
8.10 7.1872e-19 2.5524e-19 1.5363e-19 7.1872e-19 4.6924e-19
8.30 2.1006e-19 4.5208e-20 3.9567e-20 1.7857e-19 7.3093e-20
8.50 1.0553e-20 1.0553e-20 1.0585e-20 1.0553e-20 3.2064e-20
8.70 7.6933e-25 7.6933e-25 3.9934e-21 7.6933e-25 0.0000e+00
8.90 6.3612e-25 6.3612e-25 1.5606e-21 6.2234e-25 2.9150e-21
Iron
6.55 2.0076e-14 4.3579e-15 4.6002e-15 8.9445e-15 1.7332e-14
6.65 1.0708e-14 2.5729e-15 2.7444e-15 4.7453e-15 9.4494e-15
6.75 4.7949e-15 8.9370e-16 9.4942e-16 2.0718e-15 4.1329e-15
6.85 2.5608e-15 4.7747e-16 4.8946e-16 1.1053e-15 2.0607e-15
6.95 1.3541e-15 1.6494e-16 1.7640e-16 4.9651e-16 9.8859e-16
7.05 5.4684e-16 1.3947e-16 1.4277e-16 2.9046e-16 5.3660e-16
7.15 3.4501e-16 5.4846e-17 5.6698e-17 1.4362e-16 2.7911e-16
7.25 2.2919e-16 3.3237e-17 3.5190e-17 9.5254e-17 1.5320e-16
7.35 9.3776e-17 1.7256e-17 1.7585e-17 5.1983e-17 9.0199e-17
7.45 4.4824e-17 6.7859e-18 7.0274e-18 2.1827e-17 4.4183e-17
7.55 3.1772e-17 3.2806e-18 3.5274e-18 1.5352e-17 2.4856e-17
7.65 1.3672e-17 1.4549e-18 1.6190e-18 6.3086e-18 1.3534e-17
7.75 7.9038e-18 8.4219e-19 8.5326e-19 4.1903e-18 6.5355e-18
7.85 4.8619e-18 4.4056e-19 4.4546e-19 2.1417e-18 3.5520e-18
7.95 2.5898e-18 2.9273e-19 2.9924e-19 1.1319e-18 1.8093e-18
8.10 1.0928e-18 1.2142e-19 1.3512e-19 4.8093e-19 6.8341e-19
8.30 1.9928e-19 3.3365e-20 3.5448e-20 8.1524e-20 1.8182e-19
8.50 7.5703e-20 1.0430e-20 1.0916e-20 3.1102e-20 2.1113e-20
8.70 1.4712e-20 2.5145e-21 2.8635e-21 1.3300e-21 4.1937e-21
8.90 2.4809e-21 1.0360e-21 1.2659e-21 1.0434e-21 1.6930e-21

.5 K.D.E. templates and fit results

As discussed in Section III.4, the individual templates for each energy bin are shown in Figure 25, while the fit to the data for each energy bin is then shown in Figure 26. Correlation matrices for the fit results for all energy bins are shown in Table 6.

Figure 25: Distributions of the natural logarithm of the neural network mass output for each slice in energy used in the coincident analysis. Energy ranges are labeled in log{}_{10}(E{}_{0,reco}/GeV) in the titles of each figure. The y-axis represents the the number of simulated events for proton (red), helium (orange), oxygen (green) and iron (blue). The solid line represents the template probability density functions found by the adaptive KDE fitting method.
Figure 26: Best fit of the experimental data histograms by the KDE templates (derived from simulation and shown independently in Figure 25) for the coincident analysis. Energy ranges are labeled in log{}_{10}(E{}_{0,reco}/GeV) in the titles of each figure. The y-axis represents the the number of data events and the solid lines represent the weighted KDE templates for proton (red), helium (orange), oxygen (green) and iron (blue). The solid black line represents best fit distribution.
Table 6: Correlation coefficients from K.D.E. fits.
\begin{array}[]{l|l}\text{log(E/GeV)}=6.5-6.6&\text{log(E/GeV)}=6.6-6.7\\ cor=\begin{bmatrix}1.000&-0.889&0.660&-0.428\\ -0.889&1.000&-0.874&0.625\\ 0.660&-0.874&1.000&-0.878\\ -0.428&0.625&-0.878&1.000\\ \end{bmatrix}&\par cor=\begin{bmatrix}1.000&-0.786&0.544&-0.428\\ -0.786&1.000&-0.888&0.755\\ 0.544&-0.888&1.000&-0.938\\ -0.428&0.755&-0.938&1.000\\ \end{bmatrix}\\ \hline\text{log(E/GeV)}=6.7-6.8&\text{log(E/GeV)}=6.8-6.9\\ cor=\begin{bmatrix}1.000&-0.866&0.616&-0.401\\ -0.866&1.000&-0.865&0.619\\ 0.616&-0.865&1.000&-0.862\\ -0.401&0.619&-0.862&1.000\\ \end{bmatrix}&\par cor=\begin{bmatrix}1.000&-0.877&0.656&-0.455\\ -0.877&1.000&-0.887&0.670\\ 0.656&-0.887&1.000&-0.878\\ -0.455&0.670&-0.878&1.000\\ \end{bmatrix}\\ \hline\text{log(E/GeV)}=6.9-7.0&\text{log(E/GeV)}=7.0-7.1\\ cor=\begin{bmatrix}1.000&-0.778&0.466&-0.254\\ -0.778&1.000&-0.833&0.527\\ 0.466&-0.833&1.000&-0.784\\ -0.254&0.527&-0.784&1.000\\ \end{bmatrix}&\par cor=\begin{bmatrix}1.000&-0.883&0.669&-0.505\\ -0.883&1.000&-0.883&0.699\\ 0.669&-0.883&1.000&-0.891\\ -0.505&0.699&-0.891&1.000\\ \end{bmatrix}\\ \hline\text{log(E/GeV)}=7.1-7.2&\text{log(E/GeV)}=7.2-7.3\\ cor=\begin{bmatrix}1.000&-0.829&0.598&-0.378\\ -0.829&1.000&-0.884&0.613\\ 0.598&-0.884&1.000&-0.819\\ -0.378&0.613&-0.819&1.000\\ \end{bmatrix}&\par cor=\begin{bmatrix}1.000&-0.773&0.586&-0.427\\ -0.773&1.000&-0.913&0.717\\ 0.586&-0.913&1.000&-0.874\\ -0.427&0.717&-0.874&1.000\\ \end{bmatrix}\\ \hline\text{log(E/GeV)}=7.3-7.4&\text{log(E/GeV)}=7.4-7.5\\ cor=\begin{bmatrix}1.000&-0.873&0.653&-0.464\\ -0.873&1.000&-0.881&0.669\\ 0.653&-0.881&1.000&-0.870\\ -0.464&0.669&-0.870&1.000\\ \end{bmatrix}&\par cor=\begin{bmatrix}1.000&-0.889&0.660&-0.416\\ -0.889&1.000&-0.866&0.587\\ 0.660&-0.866&1.000&-0.803\\ -0.416&0.587&-0.803&1.000\\ \end{bmatrix}\\ \hline\text{log(E/GeV)}=7.5-7.6&\text{log(E/GeV)}=7.6-7.7\\ cor=\begin{bmatrix}1.000&-0.885&0.556&-0.283\\ -0.885&1.000&-0.782&0.440\\ 0.556&-0.782&1.000&-0.748\\ -0.283&0.440&-0.748&1.000\\ \end{bmatrix}&\par cor=\begin{bmatrix}1.000&-0.688&0.166&-0.014\\ -0.688&1.000&-0.647&0.279\\ 0.166&-0.647&1.000&-0.679\\ -0.014&0.279&-0.679&1.000\\ \end{bmatrix}\\ \hline\text{log(E/GeV)}=7.7-7.8&\text{log(E/GeV)}=7.8-7.9\\ cor=\begin{bmatrix}1.000&-0.789&0.438&-0.189\\ -0.789&1.000&-0.769&0.401\\ 0.438&-0.769&1.000&-0.717\\ -0.189&0.401&-0.717&1.000\\ \end{bmatrix}&\par cor=\begin{bmatrix}1.000&-0.856&0.596&-0.296\\ -0.856&1.000&-0.840&0.461\\ 0.596&-0.840&1.000&-0.712\\ -0.296&0.461&-0.712&1.000\\ \end{bmatrix}\\ \par\hline\text{log(E/GeV)}=7.9-8.0&\text{log(E/GeV)}=8.0-8.2\\ cor=\begin{bmatrix}1.000&-0.895&0.657&-0.386\\ -0.895&1.000&-0.866&0.555\\ 0.657&-0.866&1.000&-0.784\\ -0.386&0.555&-0.784&1.000\\ \end{bmatrix}&\par cor=\begin{bmatrix}1.000&-0.716&0.289&-0.035\\ -0.716&1.000&-0.764&0.318\\ 0.289&-0.764&1.000&-0.649\\ -0.035&0.318&-0.649&1.000\\ \end{bmatrix}\\ \hline\text{log(E/GeV)}=8.2-8.4&\text{log(E/GeV)}=8.4-8.6\\ cor=\begin{bmatrix}1.000&-0.013&-0.419&0.186\\ -0.013&1.000&-0.012&0.005\\ -0.419&-0.012&1.000&-0.727\\ 0.186&0.005&-0.727&1.000\\ \end{bmatrix}&\par cor=\begin{bmatrix}1.000&-0.029&-0.406&0.188\\ -0.029&1.000&-0.023&0.009\\ -0.406&-0.023&1.000&-0.790\\ 0.188&0.009&-0.790&1.000\\ \end{bmatrix}\\ \hline\text{log(E/GeV)}=8.6-8.8&\text{log(E/GeV)}=8.8-9.0\\ cor=\begin{bmatrix}1.000&-0.835&0.007&0.254\\ -0.835&1.000&-0.027&-0.469\\ 0.007&-0.027&1.000&-0.013\\ 0.254&-0.469&-0.013&1.000\\ \end{bmatrix}&\par cor=\begin{bmatrix}1.000&-0.659&0.005&0.188\\ -0.659&1.000&-0.031&-0.570\\ 0.005&-0.031&1.000&-0.024\\ 0.188&-0.570&-0.024&1.000\\ \end{bmatrix}\\ \par\end{array}
Comments 0
Request Comment
You are adding the first comment!
How to quickly get a good reply:
  • Give credit where it’s due by listing out the positive aspects of a paper before getting into which changes should be made.
  • Be specific in your critique, and provide supporting evidence with appropriate references to substantiate general statements.
  • Your comment should inspire ideas to flow and help the author improves the paper.

The better we are at sharing our knowledge with each other, the faster we move forward.
""
The feedback must be of minimum 40 characters and the title a minimum of 5 characters
   
Add comment
Cancel
Loading ...
375227
This is a comment super asjknd jkasnjk adsnkj
Upvote
Downvote
""
The feedback must be of minumum 40 characters
The feedback must be of minumum 40 characters
Submit
Cancel

You are asking your first question!
How to quickly get a good answer:
  • Keep your question short and to the point
  • Check for grammar or spelling errors.
  • Phrase it like a question
Test
Test description