RdSdSimCheck/Verification.cc

Go to the documentation of this file.

#include <evt/Event.h>
#include <evt/ShowerRecData.h>

#include <utl/ErrorLogger.h>
#include <utl/Branch.h>
#include <utl/Particle.h>
#include <utl/AugerUnits.h>

#include <fwk/CentralConfig.h>

#include <revt/REvent.h>
#include <revt/Header.h>
#include <revt/Station.h>
#include <revt/StationRecData.h>
#include <revt/StationRRecDataQuantities.h>
#include <evt/ShowerRecData.h>
#include <evt/ShowerRRecData.h>

#include <utl/Trace.h>
#include <utl/TraceAlgorithm.h>
#include <utl/ErrorLogger.h>
#include <utl/Reader.h>
#include <utl/config.h>
#include <utl/AugerUnits.h>
#include <utl/CoordinateSystemPtr.h>
#include <utl/CoordinateSystem.h>
#include <fwk/LocalCoordinateSystem.h>
#include <fwk/CoordinateSystemRegistry.h>

#include "Verification.h"
#include <tst/Validatrix.h>

#include <iostream>
#include <string>

using namespace std;
using namespace evt;
using namespace fwk;
using namespace utl;
using namespace revt;
using namespace tst::Validatrix;


namespace VerificationRadio {

  VModule::ResultFlag
  Verification::Init()
  {
    Branch topB = CentralConfig::GetInstance()->GetTopBranch("Verification");

    // get the coordinate system information
    const string csString = topB.GetChild("coordinateSystemId").Get<string>();
    fReferenceCS = fwk::CoordinateSystemRegistry::Get(csString);

    // loop through all the children and fill the list of tests to be performed
    for (Branch b = topB.GetFirstChild(); b; b = b.GetNextSibling()) {
      // instantiate and initialize variables to hold values
      unsigned long parameterId = 0;
      string parameterName;
      bool useAbsoluteTolerance = false;
      double tolerance = 0;
      // see if current branch is relevant
      if (b.GetName() == "StationQuantity") {
        b.GetChild("ParameterId").GetData(parameterId);
        b.GetChild("ParameterName").GetData(parameterName);
        b.GetChild("UseAbsoluteTolerance").GetData(useAbsoluteTolerance);
        b.GetChild("Tolerance").GetData(tolerance);
        fStationTests.push_back(TestInformation(parameterId, parameterName, useAbsoluteTolerance, tolerance));
      } else if (b.GetName() == "ShowerQuantity") {
        b.GetChild("ParameterId").GetData(parameterId);
        b.GetChild("ParameterName").GetData(parameterName);
        b.GetChild("UseAbsoluteTolerance").GetData(useAbsoluteTolerance);
        b.GetChild("Tolerance").GetData(tolerance);
        fShowerTests.push_back(TestInformation(parameterId, parameterName, useAbsoluteTolerance, tolerance));
      }
    }

    // get the name of the file hosting the reference values
    topB.GetChild("RefFileName").GetData(fRefFileName);

    // find out if we should compare or create the reference values

    // set up the stringstream for writing results
    fTestResults = new ofstream(fRefFileName + ".new");

    return fTestResults->is_open() ? eSuccess : eFailure;
  }


  VModule::ResultFlag
  Verification::Run(evt::Event& event)
  {
    ostringstream info;
    info << "Run validation for radio event E" << event.GetREvent().GetHeader().GetId();
    INFO(info);

    *fTestResults << BeLabel("Verification for Radio:");

    if (!event.HasREvent()) {
      *fTestResults << BeLabel("No radio event found! Test failed");
      ERROR("No radio event found! Test failed");
      return eFailure;
    }
    REvent& rEvent = event.GetREvent();

    if (!event.HasRecShower()) {
      *fTestResults << BeLabel("No reconstructed shower found! Test failed");
      ERROR("No reconstructed shower found! Test failed");
      return eFailure;
    }

    if (!event.GetRecShower().HasRRecShower()) {
      *fTestResults << BeLabel("No reconstructed radio shower found! Test failed");
      ERROR("No reconstructed radio shower found! Test failed");
      return eFailure;
    }
    ShowerRRecData& rRecShower=event.GetRecShower().GetRRecShower();

    // check the event id
    *fTestResults << BeLabel("Run number:")
                  << BeEqual(rEvent.GetHeader().GetRunNumber())
                  << BeLabel("Eventid:")
                  << BeEqual(rEvent.GetHeader().GetId());

    // do all the tests on shower level
    *fTestResults << BeLabel("Tests on radio shower level:");
    stringstream label;
    for (const auto& info : fShowerTests) {
      label.str("");
      label << "R" << rEvent.GetHeader().GetRunNumber() << "-E" << rEvent.GetHeader().GetId() << "-" << info.fParameterId << "->" << info.fParameterName;
      *fTestResults << BeLabel(label.str());
      if (rRecShower.HasParameter(ShowerRRecDataQuantities(info.fParameterId))) {
        if (info.fUseAbsoluteTolerance)
          *fTestResults << BeCloseAbs(rRecShower.GetParameter(ShowerRRecDataQuantities(info.fParameterId)), info.fTolerance);
        else
          *fTestResults << BeCloseRel(rRecShower.GetParameter(ShowerRRecDataQuantities(info.fParameterId)), info.fTolerance);
      } else {
        *fTestResults << BeEqual(-1);
      }
      label << "Error";
      *fTestResults << BeLabel(label.str());
      if (rRecShower.HasParameterError(ShowerRRecDataQuantities(info.fParameterId))) {
        if (info.fUseAbsoluteTolerance)
          *fTestResults << BeCloseAbs(rRecShower.GetParameterError(ShowerRRecDataQuantities(info.fParameterId)), info.fTolerance);
        else
          *fTestResults << BeCloseRel(rRecShower.GetParameterError(ShowerRRecDataQuantities(info.fParameterId)), info.fTolerance);
      } else {
        *fTestResults << BeEqual(-1);
      }
    }

    // now do all the tests on station level
    *fTestResults << BeLabel("Tests on radio station level:");
    for (REvent::StationIterator sIt = rEvent.StationsBegin(); sIt != rEvent.StationsEnd(); ++sIt) {
      Station& station = *sIt;

      if (!station.HasRecData()) {
        *fTestResults << BeLabel("No radio reconstructed data on station level found! Test failed");
        ERROR("No radio reconstructed data on station level found! Test failed");
        return eFailure;
      }

      for (const auto& info : fStationTests) {
        label.str("");
        label << "R" << rEvent.GetHeader().GetRunNumber() << "-E" << rEvent.GetHeader().GetId() << "-S" << station.GetId() << "-" << info.fParameterId << "->" << info.fParameterName;
        *fTestResults << BeLabel(label.str());
        if (station.GetRecData().HasParameter(StationRRecDataQuantities(info.fParameterId))) {
          if (info.fUseAbsoluteTolerance)
            *fTestResults << BeCloseAbs(station.GetRecData().GetParameter(StationRRecDataQuantities(info.fParameterId)), info.fTolerance);
          else
            *fTestResults << BeCloseRel(station.GetRecData().GetParameter(StationRRecDataQuantities(info.fParameterId)), info.fTolerance);
        } else {
          *fTestResults << BeEqual(-1);
        }
        label << "Error";
        *fTestResults << BeLabel(label.str());
        if (station.GetRecData().HasParameterError(StationRRecDataQuantities(info.fParameterId))) {
          if (info.fUseAbsoluteTolerance)
            *fTestResults << BeCloseAbs(station.GetRecData().GetParameterError(StationRRecDataQuantities(info.fParameterId)), info.fTolerance);
          else
            *fTestResults << BeCloseRel(station.GetRecData().GetParameterError(StationRRecDataQuantities(info.fParameterId)), info.fTolerance);
        } else {
          *fTestResults << BeEqual(-1);
        }
      }
    }

    return eSuccess;
  }


  VModule::ResultFlag
  Verification::Finish()
  {
    delete fTestResults;
    fTestResults = nullptr;

    INFO("Comparing radio reconstruction results to reference file");

    const auto comp = Compare(fRefFileName, fRefFileName + ".new", /*failOnFirst=*/false);

    return comp ? eSuccess : eFailure;
  }

}