RadioFileAERA.cc

Go to the documentation of this file.

#include <io/RadioFileAERA.h>

#include <utl/ErrorLogger.h>
#include <utl/AugerException.h>

#include <config.h>

#include <evt/Event.h>
#include <evt/Header.h>
#include <revt/REvent.h>
#include <revt/Station.h>
#include <revt/Channel.h>
#include <revt/Header.h>
#include <revt/StationTriggerData.h>
#include <revt/StationGPSData.h>
#include <revt/EventTrigger.h>
#include <revt/StationHeader.h>
#include <utl/Trace.h>
#include <utl/TimeStamp.h>
#include <utl/UTCDateTime.h>
#include <utl/LeapSeconds.h>

#include <det/Detector.h>
#include <rdet/RDetector.h>
#include <rdet/Station.h>

#include <TKey.h>
#include <TFile.h>
#include <TTree.h>

#include <cstddef>
#include <cstdlib>
#include <iostream>
#include <sstream>

#include <aevread.h>
#include <AeraEvent.h>
#include <amsg.h>
#include <llog.h>

using namespace std;
using namespace utl;
using namespace io;
using namespace evt;
using namespace revt;
using namespace aera;


RadioFileAERA::RadioFileAERA(const std::string& filename,
                             const Mode mode,
                             utl::Branch* const b)
{
  Open(filename, mode, b);
}


void
RadioFileAERA::Open(const std::string& filename, const Mode /*mode*/, utl::Branch* const /*b*/)
{
  if (fHasFileOpen)
    Close();
  setAevreadLoglevel(0); //aevread library will be silent
  setAevreadBufferedRead(true);
  if (openAeraEventFile(filename.c_str())) {
    fHasFileOpen = true;
    fCurrentEvent = 0;
  } else {
    const string msg = "Something wrong occured when opening " + filename;
    throw utl::IOFailureException(msg);
  }
}


void
RadioFileAERA::Close()
{
  if (!fHasFileOpen)
    return;
  closeAeraEventFile();
  fHasFileOpen = false;
  fCurrentEvent = -1;
}


void
RadioFileAERA::Write(const evt::Event& /*event*/)
{
  WARNING("RadioFileAERA::Write not implemented!");
}


Status
RadioFileAERA::Read(evt::Event& event)
{
  AeraEvent* const aev = fIsFirstevent ? getFirstAeraEvent() : getNextAeraEvent();
  fIsFirstevent = false;

  if (!aev) {
    fIsFirstevent = true;
    return eEOF;
  }

  if (event.HasREvent())
    WARNING("Event not cleared - has REvent");
  else
    event.MakeREvent();

  REvent& revent = event.GetREvent();
  const int event_id = getAeraEventEventId();
  const int run_id = getAeraEventRunId();
  const int file_id = getAeraFileHeaderFileId();
  //const int runningid = getAeraEventRunningEventId();
  const int eventversion = getAeraEventEventVersion();
  const int numberofstations = getAeraEventLsCount();
  revent.GetHeader().SetId(event_id); // integer ID for REvent
  revent.GetHeader().SetRunNumber(run_id);
  revent.GetHeader().SetFileId(file_id);
  //revent.GetHeader().SetEventNumber(runningid);
  revent.GetHeader().SetEventVersion(eventversion);
  revent.GetHeader().SetNumberOfStations(numberofstations);
  ostringstream id;
  id << "revt_"<< run_id<< "_" << event_id;
  event.GetHeader().SetId(id.str()); // string ID for Event
  unsigned int ts_s = getAeraEventGPSSeconds();
  unsigned int ts_ns = getAeraEventNanoSeconds();
  // fixme MM : For these 2 workaroungs a maximal runnumber has to be defined
  if (ts_s < 65536) {
    // Workaround for events with wrong timestamp
    WARNING("Event skipped! TimeStamp was not filled correctly.");
    return eSkipEvent;
  }

  if (ts_s > 1293494416) {
    //Workaround for unix timestamp bug with German data
    //Assume that everydate > 1/1/2020 is UNIX based
    // This workaround is done again around line 310
    LeapSeconds ls;
    unsigned long ts_s_gps = 0;
    ls.ConvertUnixToGPS(time_t(ts_s),ts_s_gps);
    ts_s = ts_s_gps;
  }
  utl::TimeStamp ts;
  try {
    ts.SetGPSTime(ts_s, ts_ns);
  } catch (OutOfBoundException&) {
    WARNING("Event skipped! TimeStamp was not filled correctly.");
    return eSkipEvent;
  }
  det::Detector::GetInstance().Update(ts);
  event.GetHeader().SetTime(ts);
  revent.GetHeader().SetTime(ts);

  if (!getAeraEventFirstLsData()) {
    WARNING("No Local station data found");
    return eFail;
  }

  if (revent.HasTrigger()) {
    const string msg = "Warning overwrites Event Trigger information";
    WARNING(msg);
  } else
    revent.MakeTrigger();

  if (getAeraEventIsSelfTriggered())
    revent.GetTrigger().SetSelfTrigger(true);
  if (getAeraEventIsExternalT3Trigger())
    revent.GetTrigger().SetExternalTrigger(true);
  if (getAeraEventIsExternalElTrigger())
    revent.GetTrigger().SetExternalTrigger(true);
  if (getAeraEventIsRandomTrigger())
    revent.GetTrigger().SetPeriodicTrigger(true);
  if (getAeraEventIsCalibTrigger())
    revent.GetTrigger().SetCalibrationTrigger(true);

  revent.GetTrigger().SetBitPattern(getAeraEvbTriggerFlag());
  for (int i = 0; i < getAeraEventLsCount(); ++i) { //LS station loop

    const int station_id = LSPOS(getAeraEvbLsId());

    if (!getAeraEvbChannelPresentCount())  {
      if (!getAeraEventNextLsData())
        break;
      else
        continue;
    }

    ostringstream stationname;

    if (station_id > 200 && station_id < 210)
      stationname << 'M' << (station_id - 200);
     else
      stationname << "AERA_" << station_id;
    try {
      if (revent.HasStation(stationname.str())) {
        const string msg = "Warning overwrites station information";
        WARNING(msg);
      } else
        revent.MakeStation(stationname.str());
      revt::Station& station = revent.GetStationByName(stationname.str());
      const rdet::Station& detstation =
        det::Detector::GetInstance().GetRDetector().GetStation(station); //Get the RDetector::Station from revt::Station
      if (station.HasTriggerData()) {
        const string msg = "Overwriting station trigger information";
        WARNING(msg);
      } else
        station.MakeTriggerData();
      // fixme SJ: Station also sets event level trigger because this flag is often not set correctly in the data
      if (getAeraEvbIsSelfTriggered() ||
          (getAeraEvbHardwareIsDutch() && !getAeraEvbIsRandomTrigger() && !getAeraEvbIsCalibTrigger())) {
        station.GetTriggerData().SetTriggerSource(StationTriggerData::eSelf);
        revent.GetTrigger().SetSelfTrigger(true);
      } else if (getAeraEvbIsExternalElTrigger())
        station.GetTriggerData().SetTriggerSource(StationTriggerData::eExternal);
      else if (getAeraEvbIsExternalT3Trigger()) {
        station.GetTriggerData().SetTriggerSource(StationTriggerData::eExternal);
        revent.GetTrigger().SetExternalTrigger(true);
        // WorkAround for Periodic trigger
        if ((run_id >= 100210 && run_id < 100220) && (ts_ns < 5e8+20) && (ts_ns > 5e8-20)) {
          station.GetTriggerData().SetTriggerSource(StationTriggerData::ePeriodic);
          revent.GetTrigger().SetPeriodicTrigger(true);
          revent.GetTrigger().SetExternalTrigger(false);
        }
        if (run_id >= 100220 && (ts_ns < 7 || ts_ns > 1e9-7)) {
          station.GetTriggerData().SetTriggerSource(StationTriggerData::ePeriodic);
          revent.GetTrigger().SetPeriodicTrigger(true);
          revent.GetTrigger().SetExternalTrigger(false);
        }
      } else if (getAeraEvbIsCalibTrigger())
        station.GetTriggerData().SetTriggerSource(StationTriggerData::eCalib);
      else if (getAeraEvbIsRandomTrigger())
        station.GetTriggerData().SetTriggerSource(StationTriggerData::ePeriodic);
      station.GetTriggerData().SetOffset(getAeraEvbTriggerPos());
      if (station.HasGPSData()) {
        const string msg = "Warning overwrites station GPS information";
        WARNING(msg);
      } else
        station.MakeGPSData();
      StationGPSData& gps = station.GetGPSData();

      if (station.HasStationHeader()) {
        const string msg ="Warning overwrites station header";
        WARNING(msg);
      } else
        station.MakeStationHeader();
      if (getAeraEvbHardwareIsDutch())
        station.GetStationHeader().SetHardwareType(StationHeader::eNL);
      else if (getAeraEvbHardwareIsFrench())
        station.GetStationHeader().SetHardwareType(StationHeader::eFR);
      else if (getAeraEvbHardwareIsGerman())
        station.GetStationHeader().SetHardwareType(StationHeader::eDE);
      else {
        ostringstream msg;
        msg << " Unknown Hardware type ::   getAeraEvbHardwareType() = " << getAeraEvbHardwareType() << "\n"
               " Event ID is << " << id.str() << " File is " << file_id << "\n"
               " Skip to next event";
        ERROR(msg);
        continue;
      }
      station.GetStationHeader().SetHardwareVersion(getAeraEvbHardwareVersion());

      const double samplingfreq = double(getAeraEvbSamplingFreq()) * utl::megahertz;
      ts_s = getAeraEvbGPSseconds();
      // Same Workaround as before for at the station level !!!
      if (ts_s > 1293494416) {
        // Workaround for unix timestamp bug with German data
        // Assume that everydate > 1/1/2020 is UNIX based
        unsigned long ts_s_gps = 0;
        LeapSeconds ls;
        ls.ConvertUnixToGPS(time_t(ts_s), ts_s_gps);
        ts_s = ts_s_gps;
      }
      utl::TimeStamp stationts;
      try {
        stationts.SetGPSTime(ts_s, getAeraEvbGPSnanoseconds());
        gps.SetSecond(ts_s);  // ? Why was ts_gps_s
        gps.SetCorrectedNanosecond(getAeraEvbGPSnanoseconds());  // corrected from what ?
      } catch (utl::OutOfBoundException&) {
        WARNING("Event skipped! TimeStamp was not filled correctly.");
        return eSkipEvent;
      }
      double pretrigtime = 0;
      if (getAeraEvbHardwareIsDutch()) {
        const char* const raw = (const char*)getAeraEvbInfoADCbuffer();
        const UINT16 s = (raw[5] << 8) + (raw[6] & 0xff);
        const UINT16 t = (raw[7] << 8) + (raw[8] & 0xff);
        pretrigtime = (s + t) / samplingfreq;
      } else
        pretrigtime = double(getAeraEvbTriggerPos()) / samplingfreq;
      const utl::TimeInterval pretrig = pretrigtime;
      try {
        station.SetRawTraceStartTime(stationts - pretrig);
      } catch (utl::OutOfBoundException&) {
        WARNING("Event skipped! TimeStamp was not filled correctly.");
        return eSkipEvent;
      }

      const int numberofsamples = getAeraEvbTracelength();

      UINT16 thresh = 0;
      ostringstream info;
      for (int ch_id = 1; ch_id <= detstation.GetNChannels(); ++ch_id) {
        info.str("");
        revt::Channel& chan = station.GetChannel(ch_id);
        //If Channel is not present in data, it will be set as inactive;
        if (!getAeraEvbChannelXPresent(ch_id)) {
          if (!((ch_id == 3 || ch_id == 4) && getAeraEvbHardwareIsGerman())) {
            // German data have no Channel 3/4 so this warning is not necessary
            info.str("");
            info << " Channel " << ch_id << " Not present in data\n Channel will be marked as non active";
            INFO(info);
          }
          chan.SetNotActive();
          continue;
        }
        //as dutch channels 3 & 4 are not filled with low gain data, they will be set as inactive
        if ((ch_id == 3 || ch_id == 4) && getAeraEvbHardwareIsDutch()) {
          chan.SetNotActive();
          continue;
        }
        chan.SetActive();
        revt::ChannelADCTimeSeries& thistrace = chan.GetChannelADCTimeSeries();
        if (getAeraEvbHardwareIsGerman()) {
          thistrace.SetBinning(1. / samplingfreq);
          int offset = 2048;
          if (getAevreadVersion() > 10201)
            offset = 0;

          INT16* const tmpseries = (INT16*)getAeraEvbChannelXData(ch_id); // Changed for DE data
          if (!tmpseries) {
            info.str("");
            info << " Pointer to time serie is NULL for channel " << ch_id
                 << " Channel will be marked as not active ";
            WARNING(info);
            chan.SetNotActive();
            continue;
          }
          for (int ii = 0; ii < numberofsamples; ++ii) {
            thistrace.PushBack(tmpseries[ii] + offset); //Change for DE electroncic
          }
          free(tmpseries);
        } else if (getAeraEvbHardwareIsDutch()) {
          thistrace.SetBinning(1. / samplingfreq);
          UINT16* const tmpseries =  getAeraEvbChannelXData(ch_id); // Changed for DE data
          if (!tmpseries) {
            info.str("");
            info << " Pointer to time serie is NULL for channel " << ch_id
                 << " Channel will be marked as not active ";
            WARNING(info);
            chan.SetNotActive();
            continue;
          }
          for (int ii = 0; ii < numberofsamples; ++ii) {
            thistrace.PushBack(tmpseries[ii]); //Change for DE electroncic
          }
          free(tmpseries);
        } else {
          info.str("");
          info << "Electronic readout is not implemented, Electronic type is " << getAeraEvbHardwareType()
               << " Channel is Set has inactive";
          chan.SetNotActive();
          WARNING(info);
        }
        chan.SetNyquistZone(1);
        // fixme SJ: if statement is workaround to avoid needless output of aevread for DE data, should be fixed in next version
        if (getAeraEvbHardwareIsDutch()) {
          thresh = getAeraEvbSignalThresholdX(ch_id);
          // 0xFFFF is the Error code returned by the library in case no threshold is found
          if (thresh != 0xFFFF)
            chan.SetSignalThreshold(thresh);
          thresh = getAeraEvbNoiseThresholdX(ch_id);
          if (thresh != 0xFFFF)
            chan.SetNoiseThreshold(thresh);
        }
      }
    } catch (utl::NonExistentComponentException& e) {
      ERROR(e.what());
      if (!getAeraEventNextLsData())
        break;
      else
        continue;
    }

    if (i+1 < getAeraEventLsCount())
      getAeraEventNextLsData();
  }

  ++fCurrentEvent;
  //freeEvbody();
  return eSuccess;
}


Status
RadioFileAERA::FileStatus()
{
  return fHasFileOpen ? eSuccess : eFail;
}


Status
RadioFileAERA::FindEvent(const unsigned int /*id*/)
{
  WARNING("RadioFileAERA::FindEvent not implemented!");
  return eEOF;
}


Status
RadioFileAERA::GotoPosition(const unsigned int /*position*/)
{
  WARNING("RadioFileAERA::GotoPosition not implemented!");
  return eEOF;
}


int
RadioFileAERA::GetNEvents()
{
  WARNING("RadioFileAERA::GetNEvents not implemented!");
  return 0;
}