testAiresShowerFile.cc

Go to the documentation of this file.

#include <io/AiresShowerFile.h>
#include <io/AiresUtilities.h>
#include <io/IoCodes.h>
#include <evt/Event.h>
#include <evt/ShowerSimData.h>
#include <utl/AugerUnits.h>
#include <utl/Particle.h>
#include <utl/NucleusProperties.h>
#include <utl/ShowerParticleIterator.h>
#include <fwk/CentralConfig.h>
#include <det/Detector.h>

#include <tst/Verify.h>
#include <cppunit/extensions/HelperMacros.h>

#include <iostream>

using io::AiresShowerFile;
using namespace fwk;
using namespace std;
using namespace tst;
using namespace utl;
using namespace tst;
using namespace tst;


class testAiresShowerFile : public CppUnit::TestFixture {

  CPPUNIT_TEST_SUITE(testAiresShowerFile);
  CPPUNIT_TEST(testConstruct);
  //CPPUNIT_TEST_EXCEPTION(testConstructException, utl::IOFailureException);
  CPPUNIT_TEST(testOpen);
  //CPPUNIT_TEST_EXCEPTION(testOpenOpen, utl::IOFailureException);
  //CPPUNIT_TEST_EXCEPTION(testOpenNonExisting, utl::IOFailureException);
  //CPPUNIT_TEST_EXCEPTION(testOpenForWrite, utl::IOFailureException);
  CPPUNIT_TEST(testClose);
  //CPPUNIT_TEST(testCloseClose);
  //CPPUNIT_TEST_EXCEPTION(testWrite, utl::IOFailureException);
  //CPPUNIT_TEST(testFindEvent);
  //CPPUNIT_TEST(testGotoPosition);
  CPPUNIT_TEST(testGetNEvents);
  CPPUNIT_TEST(testRead);
  CPPUNIT_TEST(testReadNucleus);
  CPPUNIT_TEST(testParticleRead);

  // split into tests for longitudinal development and energy deposit
  //
  CPPUNIT_TEST(testLongitudinalProfile);
  CPPUNIT_TEST(testdEdX);
  CPPUNIT_TEST(testSingleShowerLongitudinalProfile);
  CPPUNIT_TEST(testSingleShowerdEdX);

  /*CPPUNIT_TEST_EXCEPTION(testGetNEventsClosedException,
                           utl::IOFailureException);*/
  CPPUNIT_TEST_SUITE_END();

public:
  // shared data for tests

  AiresShowerFile* f;
  evt::Event event;

  void
  setUp()
  {
    f = new AiresShowerFile();
    event = evt::Event();
  }

  void tearDown()
  {
    delete f;
  }

  void testConstruct()
  {
    AiresShowerFile* file;
    CPPUNIT_ASSERT(file = new AiresShowerFile(TESTFILE));
    delete file;
  }

  void
  testConstructException()
  {
    Expected();
    AiresShowerFile f(TESTFILE, io::eWrite);
  }

  void testOpen() { f->Open(TESTFILE); }

  void
  testOpenNonExisting()
  {
    Expected();
    f->Open("NoSuchAiresEventFile.grdpcles");
  }

  void
  testOpenOpen()
  {
    f->Open(TESTFILE);
    Expected();
    f->Open(TESTFILE);
  }

  void
  testOpenForWrite()
  {
    Expected();
    f->Open(TESTFILE, io::eWrite);
  }

  void testClose() { f->Close(); }

  void
  testCloseClose()
  {
    f->Close();
    Expected();
    f->Close();
  }

  void
  testFindEvent()
  {
    f->Open(TESTFILE);
    CPPUNIT_ASSERT(f->FindEvent(2) != io::eFail);
    CPPUNIT_ASSERT(f->FindEvent(0) == io::eFail);
  }

  void
  testGotoPosition()
  {
    f->Open(TESTFILE);
    CPPUNIT_ASSERT(f->GotoPosition(0) != io::eFail);
    CPPUNIT_ASSERT(f->GotoPosition(2) == io::eFail);
  }

  void
  testWrite()
  {
    f->Open(TESTFILE);
    Expected();
    f->Write(event);
  }

  void
  testGetNEvents()
  {
    f->Open(TESTFILE);
    CPPUNIT_ASSERT(f->GetNEvents() == 5);
  }

  void
  testGetNEventsClosedException()
  {
    Expected();
    f->GetNEvents();
  }

  void
  testRead()
  {
    f->Open(TESTFILE);
    event = evt::Event();
    const int nEventsInFile = f->GetNEvents();
    int nEventsSeen = 0;
    int showerId = 1;
    while (f->Read(event) != io::eEOF) {
      ++nEventsSeen;

      CPPUNIT_ASSERT(event.HasSimShower());
      const evt::ShowerSimData& shower = event.GetSimShower();

      // The values used to do the checks here are set in AiresTest.inp
      //
      CPPUNIT_ASSERT(shower.GetPrimaryParticle() == Particle::eProton);
      CPPUNIT_ASSERT(Verify<CloseTo>(shower.GetEnergy(), 1.e20*eV));
      CPPUNIT_ASSERT(Verify<CloseTo>(shower.GetZenith(), 32.42*deg));
      CPPUNIT_ASSERT(Verify<CloseTo>(shower.GetMinRadiusCut(), 1.*m));   // from RLimsFile
      CPPUNIT_ASSERT(Verify<CloseTo>(shower.GetMaxRadiusCut(), 10.*km)); // from RLimsFile
      CPPUNIT_ASSERT(shower.GetShowerRunId() == string("AiresTest"));

      // fix this check.
      //CPPUNIT_ASSERT(Verify<CloseTo>(shower.GetAzimuth(),
      //                           145.59*deg-90.*deg-4.233*deg));

      CPPUNIT_ASSERT(shower.GetShowerNumber() == showerId);
      ++showerId;
      event = evt::Event();
    }

    CPPUNIT_ASSERT(nEventsInFile == nEventsSeen);
  }

  void
  testReadNucleus()
  {
    f->Open(TESTFILEIRON);
    event = evt::Event();
    f->Read(event);

    CPPUNIT_ASSERT(event.HasSimShower());
    const evt::ShowerSimData& shower = event.GetSimShower();

    CPPUNIT_ASSERT(shower.GetPrimaryParticle() ==
                   NucleusProperties::TypeCode(26, 56));
  }

  void
  testParticleRead()
  {
    CentralConfig::GetInstance(BOOTSTRAPFILE);
    f->Open(TESTFILE);

    for (int i = 0; i < f->GetNEvents(); i++) {

      event = evt::Event(); // clear event until readers do it.
      f->Read(event);

      CPPUNIT_ASSERT(event.HasSimShower());
      evt::ShowerSimData& shower = event.GetSimShower();

      shower.MakePosition(Point(0,0,0,
                                det::Detector::GetInstance().GetSiteCoordinateSystem()));

      double nParticles = 0;
      long num = 0;

      for (ShowerParticleIterator p = shower.GroundParticlesBegin();
           p != shower.GroundParticlesEnd(); ++p, ++num)
        nParticles += p->GetWeight();

      // reasonable numbers for a shower at 10^20 eV
      CPPUNIT_ASSERT(1e11 < nParticles && nParticles < 1e12);

    }
  }

  void
  testLongitudinalProfile()
  {
    CentralConfig::GetInstance(BOOTSTRAPFILE);
    f->Open(TESTFILE);

    // Just read the first event for this test.
    event = evt::Event();
    f->Read(event);

    CPPUNIT_ASSERT(event.HasSimShower());
    const evt::ShowerSimData& shower = event.GetSimShower();

    CPPUNIT_ASSERT(shower.HasLongitudinalProfile());
    const TabulatedFunction& prof = shower.GetLongitudinalProfile();

    const double zen = shower.GetZenith();

    CPPUNIT_ASSERT(Verify<CloseTo>(prof.Y(43.2583*(g/cm2)/cos(zen)),
                                      1.381535E+07)); // from AiresTest_s0001.t1291
    CPPUNIT_ASSERT(Verify<CloseTo>(prof.Y(242.240*(g/cm2)/cos(zen)),
                                      1.182769E+09)); // from AiresTest_s0001.t1291
  }

  void
  testdEdX()
  {
    CentralConfig::GetInstance(BOOTSTRAPFILE);
    f->Open(TESTFILE);

    // Just read the first event for this test.
    event = evt::Event();
    f->Read(event);
    CPPUNIT_ASSERT(event.HasSimShower());
    const evt::ShowerSimData& shower = event.GetSimShower();

    const double zen = shower.GetZenith();

    // ground depth (873.792) and number of obs levels (100) taken from AiresTest.inp
    const double auxStep = 873.792/100./cos(zen);

    CPPUNIT_ASSERT(shower.HasdEdX());
    const TabulatedFunction& dEdX = shower.GetdEdX();

    CPPUNIT_ASSERT(Verify<CloseTo>(dEdX.Y(43.2583*(g/cm2)/cos(zen))*auxStep/(GeV/(g/cm2)),
                                      531851.)); // from AiresTest_s0001.t7991
    CPPUNIT_ASSERT(Verify<CloseTo>(dEdX.Y(795.930*(g/cm2)/cos(zen))*auxStep/(GeV/(g/cm2)),
                                      1.494458E+09)); // from AiresTest_s0001.t7991
  }

  void
  testSingleShowerLongitudinalProfile()
  {
    f->Open(TESTFILESINGLESHOWER);
    event = evt::Event();
    f->Read(event);

    CPPUNIT_ASSERT(event.HasSimShower());
    const evt::ShowerSimData& shower = event.GetSimShower();

    CPPUNIT_ASSERT(shower.HasLongitudinalProfile());
    const TabulatedFunction& prof = shower.GetLongitudinalProfile();

    const double zen = shower.GetZenith();

    CPPUNIT_ASSERT(Verify<CloseTo>(prof.Y(43.2583*(g/cm2)/cos(zen)),
                                      182081.));      // from AiresTestSingleShower.t1291
    CPPUNIT_ASSERT(Verify<CloseTo>(prof.Y(242.240*(g/cm2)/cos(zen)),
                                      5.742057E+08)); // from AiresTestSingleShower.t1291
  }

  void
  testSingleShowerdEdX()
  {
    f->Open(TESTFILESINGLESHOWER);
    event = evt::Event();
    f->Read(event);

    CPPUNIT_ASSERT(event.HasSimShower());
    const evt::ShowerSimData& shower = event.GetSimShower();

    const double zen = shower.GetZenith();

    // ground depth (873.792) and number of obs levels (100) taken from AiresTest.inp
    const double auxStep = 873.792/100./cos(zen);

    CPPUNIT_ASSERT(shower.HasdEdX());
    const TabulatedFunction& dEdX = shower.GetdEdX();

    CPPUNIT_ASSERT(Verify<CloseTo>(dEdX.Y(43.2583*(g/cm2)/cos(zen))*auxStep/(GeV/(g/cm2)),
                                      263848.)); // from AiresTest_s0001.t7991
    CPPUNIT_ASSERT(Verify<CloseTo>(dEdX.Y(795.930*(g/cm2)/cos(zen))*auxStep/(GeV/(g/cm2)),
                                      1.598191E+09)); // from AiresTest_s0001.t7991
  }

};


CPPUNIT_TEST_SUITE_REGISTRATION(testAiresShowerFile);


// Configure (x)emacs for this file ...
// Local Variables:
// mode:c++
// compile-command: "make -C .. -k testAires && (cd ..; testAires)"
// End: