#include <cmath>
#include <cstdio>
#include <string>
#include "radixbug/bug.hh"
#include "radixio/eafstream.hh"
#include "radixio/hysplitcdump.hh"
namespace radix
{
template <typename data_type>
HysplitCDumpStream<data_type>::HysplitCDumpStream(data_type container)
{
  mData = container;
}

template <typename data_type>
bool HysplitCDumpStream<data_type>::read_from(const std::string &file)
{
  bool result = false;
  eafstream fstr(file.c_str(), std::ios::in | std::ios::binary);
  if (fstr.is_open() == false)
  {
    return result;
  }
  // big endian
  fstr.setReverseBytes(true);
  // read fortran header
  fstr.readHeader();
  std::string id = fstr.readString(4);
  mData->setId(id);
  int packing, year, month, day, hour, forecastHour, numLocations;
  fstr >> year >> month >> day >> hour >> forecastHour >> numLocations >>
      packing;
  fstr.readFooter();
  radix_ensure(fstr.header() == fstr.footer());
  mData->setDateTime(year, month, day, hour);
  mData->setForecastHour(forecastHour);

  float lat, lon, z;
  int minutes;
  // for the number of locations
  for (int i = 0; i < numLocations; ++i)
  {
    fstr.readHeader();
    fstr >> year >> month >> day >> hour >> lat >> lon >> z >> minutes;
    mData->addLocation(year, month, day, hour, lat, lon, z, minutes);
    fstr.readFooter();
    radix_ensure(fstr.header() == fstr.footer());
  }

  int numLat, numLon;
  float deltaLat, deltaLon, cLat, cLon;
  fstr.readHeader();
  fstr >> numLat >> numLon >> deltaLat >> deltaLon >> cLat >> cLon;
  fstr.readFooter();
  radix_ensure(fstr.header() == fstr.footer());
  mData->setNumLat(numLat);
  mData->setNumLon(numLon);
  mData->setDeltaLat(deltaLat);
  mData->setDeltaLon(deltaLon);
  mData->setCornerLat(cLat);
  mData->setCornerLon(cLon);

  fstr.readHeader();
  int numLevels;
  fstr >> numLevels;
  for (int i = 0; i < numLevels; ++i)
  {
    int level;
    fstr >> level;
    mData->addLevel(level);
  }
  fstr.readFooter();
  radix_ensure(fstr.header() == fstr.footer());

  fstr.readHeader();
  int numPollutants;
  fstr >> numPollutants;
  for (int i = 0; i < numPollutants; ++i)
  {
    std::string pol;
    pol = fstr.readString(4);
    mData->addPollutant(pol);
  }
  fstr.readFooter();
  radix_ensure(fstr.header() == fstr.footer());

  // infinite loop over observation times
  int timeIndex = 0;
  for (; fstr.good();)
  {
    fstr.readHeader();
    fstr >> year >> month >> day >> hour >> minutes >> forecastHour;
    fstr.readFooter();
    radix_ensure(fstr.header() == fstr.footer());
    mData->addStartTime(year, month, day, hour, minutes, forecastHour);

    fstr.readHeader();
    fstr >> year >> month >> day >> hour >> minutes >> forecastHour;
    fstr.readFooter();
    radix_ensure(fstr.header() == fstr.footer());
    mData->addEndTime(year, month, day, hour, minutes, forecastHour);

    for (int i = 0; i < numPollutants; ++i)
    {
      for (int j = 0; j < numLevels; ++j)
      {
        fstr.readHeader();
        std::string pollutant;
        int z;
        int nxyp;
        if (packing)
        {
          pollutant = fstr.readString(4);
          fstr >> z >> nxyp;
          mData->setNumNonZeroPoints(timeIndex, i, j, nxyp);
          for (int np = 0; np < nxyp; ++np)
          {
            short xi = fstr.readShort();
            short yi = fstr.readShort();
            float c  = fstr.readFloat();
            mData->setPoint(timeIndex, i, j, np, xi, yi, c);
          }
        }  // if packing is enabled
        fstr.readFooter();
        radix_ensure(fstr.header() == fstr.footer());
        fstr.peek();  // peek ahead to check file integrity
      }               // loop over num levels
    }                 // loop num pollutants
    timeIndex++;
  }  // infinite loop over output times
  result = true;
  return result;
}  // HysplitCDump::read_from

template <typename data_type>
bool HysplitCDumpStream<data_type>::write_to(const std::string &file) const
{
  bool result = false;
  eafstream fstr(file.c_str(), std::ios::out | std::ios::binary);
  if (fstr.is_open() == false)
  {
    return result;
  }
  // big endian
  fstr.setReverseBytes(true);
  int record_length;

  std::string id = mData->id();
  {
    // force string to be no more than 4 characters
    char cid[5];
    snprintf(cid, sizeof(cid), "%4s", id.c_str());
    id = cid;
  }
  record_length = id.size() + sizeof(int) * 7;
  int year, month, day, hour, forecastHour, minutes, numLocations, packing = 1;
  mData->dateTime(year, month, day, hour);
  //
  // Adjust the date from YYYY to YY format.
  // HYSPLIT utils always assume YY format
  if (year > 99)
  {
    year = year % 100;
  }
  forecastHour = mData->forecastHour();
  numLocations = mData->numLocations();

  radix_tagged_line("record_length=" << record_length);
  fstr.writeHeader(record_length);
  fstr.writeString(id);
  fstr << year << month << day << hour << forecastHour << numLocations
       << packing;
  radix_tagged_line("bytes_written=" << fstr.bytesWritten());
  fstr.writeFooter(record_length);

  float lat, lon, z;
  for (int i = 0; i < numLocations; ++i)
  {
    mData->location(i, year, month, day, hour, lat, lon, z, minutes);  //
    // Adjust the date from YYYY to YY format.
    // HYSPLIT utils always assume YY format
    if (year > 99)
    {
      year = year % 100;
    }
    record_length = sizeof(int) * 5 + sizeof(float) * 3;
    radix_tagged_line("record_length=" << record_length);
    fstr.writeHeader(record_length);
    fstr << year << month << day << hour << lat << lon << z << minutes;
    radix_tagged_line("bytes_written=" << fstr.bytesWritten());
    fstr.writeFooter(record_length);
  }
  int numLat, numLon;
  float deltaLat, deltaLon, cLat, cLon;
  numLat   = mData->numLat();
  numLon   = mData->numLon();
  deltaLat = mData->deltaLat();
  deltaLon = mData->deltaLon();
  cLat     = mData->cornerLat();
  cLon     = mData->cornerLon();

  record_length = sizeof(int) * 2 + sizeof(float) * 4;
  radix_tagged_line("record_length=" << record_length);
  fstr.writeHeader(record_length);
  fstr << numLat << numLon << deltaLat << deltaLon << cLat << cLon;
  radix_tagged_line("bytes_written=" << fstr.bytesWritten());
  fstr.writeFooter(record_length);

  int numLevels = mData->numLevels();
  record_length = int(sizeof(int)) + int(sizeof(int)) * numLevels;
  radix_tagged_line("record_length=" << record_length);
  fstr.writeHeader(record_length);
  fstr << numLevels;
  for (int i = 0; i < numLevels; ++i)
  {
    fstr << mData->level(i);
  }
  radix_tagged_line("bytes_written=" << fstr.bytesWritten());
  fstr.writeFooter(record_length);

  int numPollutants = mData->numPollutants();
  record_length     = int(sizeof(int)) + numPollutants * 4;
  radix_tagged_line("record_length=" << record_length);
  fstr.writeHeader(record_length);
  fstr << numPollutants;
  std::vector<std::string> local_pols;
  for (int i = 0; i < numPollutants; ++i)
  {
    std::string pol = mData->pollutant(i);
    radix_tagged_line("pollutant=" << pol);
    // force string to be no more than 4 characters
    char icid[5];
    snprintf(icid, sizeof(icid), "%4s", pol.c_str());
    std::string tmp = icid;
    fstr.writeString(tmp);
    local_pols.push_back(tmp);
  }
  radix_tagged_line("bytes_written=" << fstr.bytesWritten());
  fstr.writeFooter(record_length);

  int numStartTimes = mData->numStartTimes();
  int numEndTimes   = mData->numEndTimes();
  radix_insist(
      numStartTimes == numEndTimes,
      "The number of start times did not match the number of end times.");
  for (int timeIndex = 0; timeIndex < numStartTimes; ++timeIndex)
  {
    mData->startTime(timeIndex, year, month, day, hour, minutes,
                     forecastHour);  //
    // Adjust the date from YYYY to YY format.
    // HYSPLIT utils always assume YY format
    if (year > 99)
    {
      year = year % 100;
    }
    record_length = sizeof(int) * 6;
    radix_tagged_line("record_length=" << record_length);
    fstr.writeHeader(record_length);
    fstr << year << month << day << hour << minutes << forecastHour;
    radix_tagged_line("bytes_written=" << fstr.bytesWritten());
    fstr.writeFooter(record_length);

    mData->endTime(timeIndex, year, month, day, hour, minutes,
                   forecastHour);  //
    // Adjust the date from YYYY to YY format.
    // HYSPLIT utils always assume YY format
    if (year > 99)
    {
      year = year % 100;
    }
    record_length = sizeof(int) * 6;
    radix_tagged_line("record_length=" << record_length);
    fstr.writeHeader(record_length);
    fstr << year << month << day << hour << minutes << forecastHour;
    radix_tagged_line("bytes_written=" << fstr.bytesWritten());
    fstr.writeFooter(record_length);

    for (int i = 0; i < numPollutants; ++i)
    {
      for (int j = 0; j < numLevels; ++j)
      {
        int level     = mData->level(j);
        int nxyp      = mData->numNonZeroPoints(timeIndex, i, j);
        record_length = 4 + int(sizeof(int)) * 2 +
                        int(sizeof(short)) * 2 * nxyp +
                        int(sizeof(float)) * nxyp;
        radix_tagged_line("record_length="
                          << record_length << " pol=" << local_pols[size_t(i)]
                          << " z=" << level << " nxyp=" << nxyp);
        fstr.writeHeader(record_length);
        fstr.writeString(local_pols[size_t(i)]);
        fstr << level << nxyp;
        for (int np = 0; np < nxyp; ++np)
        {
          short xi, yi;
          float c;
          mData->point(timeIndex, i, j, np, xi, yi, c);
          fstr << xi << yi << c;
        }
        radix_tagged_line("bytes_written=" << fstr.bytesWritten());
        fstr.writeFooter(record_length);
      }
    }
  }
  fstr.close();
  result = true;
  return result;
}  // namespace radix

}  // namespace radix