gpxtools/gpxwind.cpp

#include <cmath>
#include <cstring>
#include <fstream>
#include <iomanip>
#include <iostream>
#include <limits>
#include <vector>

#include "Algorithm.h"
#include "Arguments.h"
#include "XMLParser.h"


namespace gpxtools
{
  class GPXWind : public XMLParserHandler
  {
  public:
    // -- Constructor -----------------------------------------------------------
    GPXWind() :
      _arguments("gpxwind [OPTION].. [N|NE|E|SE|S|SW|W|NW] [FILE]\nShow the impact of the wind direction on the track or route in the GPX-file.\n", "gpxwind v0.1", "Show the wind impact on a track or route"),
      _modeArgument(_arguments,  true,  'm', "mode",  "summary|full|plot", "mode for showing the impact", "summary"),
      _xmlParser(this)
    {
    }

    // -- Deconstructor----------------------------------------------------------
    virtual ~GPXWind()
    {
    }

    // -- Parse arguments -------------------------------------------------------
    bool parseArguments(int argc, char *argv[])
    {
      std::vector<std::string> filenames;

      if (!_arguments.parse(argc,argv, filenames)) return false;

      if (!checkArguments()) return false;

      if (filenames.empty())
      {
        std::cerr << "Missing wind direction" << std::endl;
        _arguments.printHelp();
        return false;
      }

      if (!parseWindDirection(filenames.front())) return false;

      filenames.erase(filenames.begin());

      if (filenames.empty())
      {
        return _xmlParser.parse(std::cin);
      }
      else
      {
        if (!parseFile(filenames.front())) return false;
      }

      return true;
    }

  private:
    // -- Check arguments ---------------------------------------------------------
    bool checkArguments()
    {
      if (_modeArgument.value() == "summary")
      {
        _mode = SUMMARY;
      }
      else if (_modeArgument.value() == "full")
      {
        _mode = FULL;
      }
      else if (_modeArgument.value() == "plot")
      {
        _mode = PLOT;
      }
      else
      {
        std::cerr << "Invalid value for --" << _modeArgument.longOption() << ": " << _modeArgument.value() << std::endl;
        return false;
      }

      return true;
    }

    bool parseWindDirection(const std::string &direction)
    {
      if (direction == "N")
      {
        _wind = 180.0;
      }
      else if (direction == "NE")
      {
        _wind = 225.0;
      }
      else if (direction == "E")
      {
        _wind = 270.0;
      }
      else if (direction == "SE")
      {
        _wind = 315.0;
      }
      else if (direction == "S")
      {
        _wind = 0.0;
      }
      else if (direction == "SW")
      {
        _wind = 45.0;
      }
      else if (direction == "W")
      {
        _wind = 90.0;
      }
      else if (direction == "NW")
      {
        _wind = 135.0;
      }
      else
      {
        std::cerr << "Invalid value for the wind direction:" << direction << std::endl;
        return false;
      }

      return true;
    }

    // -- Parse a file ----------------------------------------------------------
    bool parseFile(const std::string &filename)
    {
      bool ok =false;

      std::ifstream file(filename);

      if (file.is_open())
      {
        ok = _xmlParser.parse(file);

        file.close();
      }
      else
      {
        std::cerr << "Unable to open: " << filename << std::endl;
      }

      return ok;
    }

  public:
    // -- Output the result -----------------------------------------------------
    void report()
    {
      switch(_mode)
      {
        case SUMMARY:
          reportSummary();
          break;
        case FULL:
          reportFull();
          break;
        case PLOT:
          reportPlot();
          break;
      }
    }

    // -- Callbacks -------------------------------------------------------------
    virtual void startElement(const std::string &path, const std::string &, const Attributes &attributes)
    {
      if (path == "/gpx/rte")
      {
        _type = "Route";
        _name.clear();
        _number = 0;

        _segment._points.clear();
      }
      else if (path == "/gpx/rte/rtept")
      {
        _point._lat = getDoubleAttribute(attributes, "lat");
        _point._lon = getDoubleAttribute(attributes, "lon");
      }
      else if (path == "/gpx/trk")
      {
        _type = "Track";
        _name.clear();
        _number = 0;

        _segment._points.clear();
      }
      else if (path == "/gpx/trk/trkseg")
      {
        _number++;
        _segment._number = _number;

        _segment._points.clear();
      }
      else if (path == "/gpx/trk/trkseg/trkpt")
      {
        _point._lat = getDoubleAttribute(attributes, "lat");
        _point._lon = getDoubleAttribute(attributes, "lon");
      }
    }

    virtual void text(const std::string &path, const std::string &text)
    {
      if (path == "/gpx/rte/name")
      {
        _name = text;
      }
      else if (path == "/gpx/rte/desc")
      {
        if (_name.empty()) _name = text;
      }
      else if (path == "/gpx/trk/name")
      {
        _name = text;
      }
      else if (path == "/gpx/trk/desc")
      {
        if (_name.empty()) _name = text;
      }
    }

    virtual void endElement(const std::string &path, const std::string &)
    {
      if (path == "/gpx/rte/rtept")
      {
        _segment._points.push_back(_point);
      }
      else if (path == "/gpx/rte")
      {
        _segment._name = _name;
        _segment._type = _type;

        _segments.push_back(_segment);
      }
      else if (path == "/gpx/trk")
      {
      }
      else if (path == "/gpx/trk/trkseg")
      {
        _segment._name = _name;
        _segment._type = _type;

        _segments.push_back(_segment);
      }
      else if (path == "/gpx/trk/trkseg/trkpt")
      {
        _segment._points.push_back(_point);
      }
    }

  // -- Privates ----------------------------------------------------------------
  private:
    struct LatLon
    {
      double _lat = NAN;
      double _lon = NAN;
    };

    struct Bounds
    {
      LatLon _min = {std::numeric_limits<double>::max(), std::numeric_limits<double>::max()};
      LatLon _max = {std::numeric_limits<double>::min(), std::numeric_limits<double>::min()};

      void add(const LatLon &point)
      {
        _min._lat = std::min(point._lat, _min._lat);
        _min._lon = std::min(point._lon, _min._lon);
        _max._lat = std::max(point._lat, _max._lat);
        _max._lon = std::max(point._lon, _max._lon);
      }
    };

    struct Line
    {
      double _distance = NAN;
      double _bearing = NAN;
    };

    struct Segment
    {
      std::string _type;
      std::string _name;
      unsigned _number = 0;

      std::vector<LatLon> _points;
    };

    double getDouble(const std::string &value)
    {
      try
      {
        if (value.empty())
        {
          return NAN;
        }
        else
        {
          return std::stod(value);
        }
      }
      catch (...)
      {
        return NAN;
      }
    }

    double getDoubleAttribute(const Attributes &atts, const std::string &key)
    {
      auto iter = atts.find(key);

      return iter != atts.end() ? getDouble(iter->second) : NAN;
    }

    bool getLine(LatLon &prev, const LatLon &point, Line &line)
    {
      bool ok = false;

      if ((!std::isnan(prev._lat)) && (!std::isnan(prev._lon)) && (!std::isnan(point._lat)) && (!std::isnan(point._lon)))
      {
        line._distance = gpx::calcDistance(prev._lat, prev._lon, point._lat, point._lon);
        line._bearing = gpx::calcBearingInDeg(prev._lat, prev._lon, point._lat, point._lon);

        ok = true;
      }

      prev._lat = point._lat;
      prev._lon = point._lon;

      return ok;
    }

    void reportHeader(const Segment &segment)
    {
      std::cout << segment._type;
      if (!segment._name.empty()) std::cout << ": " << segment._name;
      if (segment._number > 0) std::cout << " (" << segment._number << ")";
      std::cout << std::endl;
    }

    void reportSummary()
    {
      for (const auto &segment : _segments)
      {
        double headwind = 0.0;
        double tailwind = 0.0;

        LatLon prev;
        Line line;

        for (const auto &point : segment._points)
        {
          if (!getLine(prev, point, line)) continue;

          double wind = std::abs(line._bearing - _wind);

          if (wind > 180.0) wind = 360.0 - wind; // Range 0..180.0

          if (wind > 90.0) // 90.0..180.0
          {
            headwind += line._distance;
          }
          else // 0..90
          {
            tailwind += line._distance;
          }
        }

        reportHeader(segment);
        std::cout << "  Headwind: " << std::fixed << std::setw(8) << std::setprecision(3) << headwind / 1000.0 << " km" << std::endl;
        std::cout << "  Tailwind: " << std::fixed << std::setw(8) << std::setprecision(3) << tailwind / 1000.0 << " km" << std::endl;
      }
    }

    void reportFull()
    {
      for (const auto &segment : _segments)
      {
        double headwind = 0.0;
        double headcrosswind = 0.0;
        double crosswind = 0.0;
        double tailcrosswind = 0.0;
        double tailwind = 0.0;

        LatLon prev;
        Line line;

        for (const auto &point : segment._points)
        {
          if (!getLine(prev, point, line)) continue;

          double wind = std::abs(line._bearing - _wind);

          if (wind > 180.0) wind = 360.0 - wind; // Range 0..180.0

          if (wind >= 170.0)
          {
            headwind += line._distance;
          }
          else if (wind > 100.0)
          {
            headcrosswind += line._distance;
          }
          else if (wind >= 80.0)
          {
            crosswind += line._distance;
          }
          else if (wind > 10.0)
          {
            tailcrosswind += line._distance;
          }
          else
          {
            tailwind += line._distance;
          }
        }

        reportHeader(segment);
        std::cout << "  Headwind      : " << std::fixed << std::setw(8) << std::setprecision(3) << headwind / 1000.0      << " km" << std::endl;
        std::cout << "  Headcrosswind : " << std::fixed << std::setw(8) << std::setprecision(3) << headcrosswind / 1000.0 << " km" << std::endl;
        std::cout << "  Crosswind     : " << std::fixed << std::setw(8) << std::setprecision(3) << crosswind / 1000.0     << " km" << std::endl;
        std::cout << "  Tailcrosswind : " << std::fixed << std::setw(8) << std::setprecision(3) << tailcrosswind / 1000.0 << " km" << std::endl;
        std::cout << "  Tailwind      : " << std::fixed << std::setw(8) << std::setprecision(3) << tailwind / 1000.0      << " km" << std::endl;
      }
    }

    void reportPlot()
    {
      Bounds bounds = reportPlotData(std::cout);

      std::cout << "set format x \"%D %E\" geographic" << std::endl;
      std::cout << "set format y \"%D %N\" geographic" << std::endl;
      std::cout << "unset key" << std::endl;
      std::cout << "unset parametric" << std::endl;
      std::cout << "unset border" << std::endl;
      std::cout << "unset xtics" << std::endl;
      std::cout << "unset ytics" << std::endl;
      std::cout << "unset colorbox" << std::endl;
      std::cout << "set style data lines" << std::endl;
      std::cout << "set xrange [ " << std::fixed << std::setprecision(6) << bounds._min._lon << " : " << bounds._max._lon << " ] noreverse nowriteback" << std::endl;
      std::cout << "set yrange [ " << std::fixed << std::setprecision(6) << bounds._min._lat << " : " << bounds._max._lat << " ] noreverse nowriteback" << std::endl;
      std::cout << "set cbrange [0.0:180.0]" << std::endl;
      std::cout << "set palette defined (0.0 \"green\", 90.0 \"yellow\", 180.0 \"red\")" << std::endl;
      std::cout << "set output \"gpxwind.png\"" << std::endl;
      reportPlotSize(std::cout, bounds);

      std::cout << "plot $DATA using 1:2:3 with lines linewidth 3 linecolor palette" << std::endl;
    }

    Bounds reportPlotData(std::ostream &str)
    {
      str << "$DATA << EOD" << std::endl;

      Bounds bounds;

      for (const auto &segment : _segments)
      {
        LatLon prev;
        Line line;

        for (const auto &point : segment._points)
        {
          if (!getLine(prev, point, line)) continue;

          double wind = std::abs(line._bearing - _wind);

          if (wind > 180.0) wind = 360.0 - wind;

          str << std::fixed << std::setprecision(6) << point._lon << " " << point._lat << " " << std::setprecision(1) << wind << std::endl;

          bounds.add(point);
        }
      }
      str << "EOD" << std::endl;

      return bounds;
    }

    void reportPlotSize(std::ostream &str, const Bounds &bounds)
    {
      const int pngSize = 1024;

      int xdistance = std::max(std::abs(static_cast<int>(gpx::calcDistance(bounds._min._lat, bounds._min._lon, bounds._min._lat, bounds._max._lon))), 1);
      int ydistance = std::max(std::abs(static_cast<int>(gpx::calcDistance(bounds._min._lat, bounds._min._lon, bounds._max._lat, bounds._min._lon))), 1);

      if (xdistance > ydistance)
      {
        ydistance = pngSize * ydistance / xdistance;
        xdistance = pngSize;
      }
      else
      {
        xdistance = pngSize * xdistance / ydistance;
        ydistance = pngSize;
      }

      str << "set term png size " << xdistance << ',' << ydistance << std::endl;
    }

    // -- Members ---------------------------------------------------------------
    arg::Arguments _arguments;
    arg::Argument _modeArgument;

    // Arguments
    enum Mode
    {
      SUMMARY, FULL, PLOT
    };
    Mode _mode = SUMMARY;

    double _wind = 0.0;

    // State of the parsing
    XMLParser _xmlParser;
    std::string _path;
    std::string _type;
    std::string _name;
    unsigned _number = 0;
    LatLon _point;
    Segment _segment;
    std::vector<Segment> _segments;
  };
};

// -- Main program ------------------------------------------------------------

int main(int argc, char *argv[])
{
  gpxtools::GPXWind gpxWind;

  if (gpxWind.parseArguments(argc, argv))
  {
    gpxWind.report();

    return EXIT_SUCCESS;
  }

  return EXIT_FAILURE;
}