Predict Example Source Files

Table of Contents

• Predict Example Source Files

  • Main.cpp

  • Predict.h

  • Predict.cpp

  • PredictData.h

  • BuildFeatures.h

  • BuildFeatures.cpp

Main.cpp

/*
 * Copyright (c) 2014-2023 National Technology and Engineering
 * Solutions of Sandia, LLC. Under the terms of Contract DE-NA0003525
 * with National Technology and Engineering Solutions of Sandia, LLC,
 * the U.S. Government retains certain rights in this software.
 *
 * Redistribution and use in source and binary forms, with or without
 * modification, are permitted provided that the following conditions
 * are met:
 *
 * 1. Redistributions of source code must retain the above copyright
 * notice, this list of conditions and the following disclaimer.
 *
 * 2. Redistributions in binary form must reproduce the above copyright
 * notice, this list of conditions and the following disclaimer in the
 * documentation and/or other materials provided with the distribution.
 *
 * THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
 * "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
 * LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR
 * A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT
 * HOLDER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL,
 * SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT
 * LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE,
 * DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY
 * THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT
 * (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE
 * OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
 */

#include "Predict.h"

#include <tracktable/CommandLineFactories/AssemblerFromCommandLine.h>
#include <tracktable/CommandLineFactories/PointReaderFromCommandLine.h>
#include <tracktable/Domain/Terrestrial.h>

#include <boost/timer/timer.hpp>

#include <string>
#include <vector>

using TrajectoryT = tracktable::domain::terrestrial::trajectory_type;
using PointT = typename TrajectoryT::point_type;
using PointReaderT = tracktable::PointReader<PointT>;
using PointReaderIteratorT = typename PointReaderT::iterator;
using AssemblerT = tracktable::AssembleTrajectories<TrajectoryT, PointReaderIteratorT>;

bool has_tail_number(const TrajectoryT &_trajectory);
bool has_destination(const TrajectoryT &_trajectory);

static constexpr auto helpmsg = R"(
--------------------------------------------------------------------------------
This example demonstrates using feature vectors to measure similarities between
trajectories via an Rtree

The predict example demonstrates:
    - Using command line factories to read points and assemble trajectories
    - Using boost program options to take parameters from command lines(in addition to the factories)
    - Conditioning trajectories based on length and objectid
    - Using boost rtree to locate similar trajectories based on cartesian distance in feature space

Typical use: '--string-field=dest x' is required

    ./predict --input=/data/SampleASDI.csv --delimiter=, --string-field=dest 30 --num-samples=10

--------------------------------------------------------------------------------)";
int main(int _argc, char *_argv[]) {
    auto numSamples = 10u;

    // Set log level to reduce unecessary output
    tracktable::set_log_level(tracktable::log::info);
    // Create a basic command line option with boost
    bpo::options_description commandLineOptions("Options");
    // clang-format off
    commandLineOptions.add_options()
      ("help", "Print help")
      ("num-samples", bpo::value(&numSamples)->default_value(10),
        "Number of samples")
    ;
    // clang-format on

    // Create command line factories
    tracktable::PointReaderFromCommandLine<PointT> readerFactory;
    tracktable::AssemblerFromCommandLine<TrajectoryT> assemblerFactory;
    // Add options from the factories
    readerFactory.addOptions(commandLineOptions);
    assemblerFactory.addOptions(commandLineOptions);

    /** Boost program options using a variable map to tie everything together.
     * one parse will have a single variable map. We need to let the factories know
     * about this variable map so they can pull information out of it */
    auto vm = std::make_shared<boost::program_options::variables_map>();
    readerFactory.setVariables(vm);
    assemblerFactory.setVariables(vm);

    // Parse the command lines, don't forget the 'notify' after
    try {
        // We use this try/catch to automatically display help when an unknown option is used
        boost::program_options::store(
            boost::program_options::command_line_parser(_argc, _argv).options(commandLineOptions).run(), *vm);
        boost::program_options::notify(*vm);
    } catch (boost::program_options::error e) {
        std::cerr << e.what();
        std::cerr << helpmsg << "\n\n";
        std::cerr << commandLineOptions << std::endl;
        return 1;
    }
    /** Parsing will give an error of an incorrect option is used, but it won't
     * display the help unless we tell it too */
    if (vm->count("help") != 0) {
        std::cerr << helpmsg << "\n\n";
        std::cerr << commandLineOptions << std::endl;
        return 1;
    }

    // Create Point Reader and assembler
    auto pointReader = readerFactory.createPointReader();
    auto assembler = assemblerFactory.createAssembler(pointReader);

    std::vector<std::shared_ptr<TrajectoryT>> trajectories = {};
    // This block exists for easy timing of trajectory assembling using the boost auto timer
    // Note that all feedback to the user is done on std::cerr, this allows us to only
    // put desired results into std::cout, this make processing output easier.
    {
        std::cerr << "Assemble Trajectories" << std::endl;
        boost::timer::auto_cpu_timer timer3(std::cerr);
        auto count = 0u;
        std::cerr << std::right;
        for (auto tIter = assembler->begin(); tIter != assembler->end(); ++tIter) {
            if (has_tail_number(*tIter)) {  // Skip tail number flights
                continue;
            }
            if (!has_destination(*tIter)) {  // Skip flights without destination or errors
                continue;
            }
            std::cerr << "\b\b\b\b\b\b\b\b\b\b" << std::setw(10)  // Using backspaces for in place counter
                      << count++;
            trajectories.push_back(std::make_shared<TrajectoryT>(*tIter));
        }
        std::cerr << std::left << "\nStarting with " << trajectories.size() << " trajectories" << std::endl;
    }

    // This routine does a prediction based on destination airport
    Predict(trajectories, numSamples);

    return 0;
}

bool has_tail_number(const TrajectoryT &_trajectory) {
    auto s = _trajectory.object_id();
    return !((s[0] != 'N') || (s[1] < '0') || (s[1] > '9'));
}

bool has_destination(const TrajectoryT &_trajectory) {
    return !_trajectory.front().string_property("dest").empty() &&
           (_trajectory.front().string_property("dest") == _trajectory.back().string_property("dest"));
}

---

Predict.h

/*
 * Copyright (c) 2014-2023 National Technology and Engineering
 * Solutions of Sandia, LLC. Under the terms of Contract DE-NA0003525
 * with National Technology and Engineering Solutions of Sandia, LLC,
 * the U.S. Government retains certain rights in this software.
 *
 * Redistribution and use in source and binary forms, with or without
 * modification, are permitted provided that the following conditions
 * are met:
 *
 * 1. Redistributions of source code must retain the above copyright
 * notice, this list of conditions and the following disclaimer.
 *
 * 2. Redistributions in binary form must reproduce the above copyright
 * notice, this list of conditions and the following disclaimer in the
 * documentation and/or other materials provided with the distribution.
 *
 * THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
 * "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
 * LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR
 * A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT
 * HOLDER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL,
 * SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT
 * LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE,
 * DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY
 * THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT
 * (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE
 * OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
 */

#ifndef Predict_h
#define Predict_h

#include <tracktable/Domain/Terrestrial.h>

#include <vector>

using TrajectoryT = tracktable::domain::terrestrial::trajectory_type;
using TrajectoryVectorT = std::vector<std::shared_ptr<TrajectoryT>>;

void Predict(const TrajectoryVectorT &_trajectories, const size_t _numSamples);

#endif

---

Predict.cpp

/*
 * Copyright (c) 2014-2023 National Technology and Engineering
 * Solutions of Sandia, LLC. Under the terms of Contract DE-NA0003525
 * with National Technology and Engineering Solutions of Sandia, LLC,
 * the U.S. Government retains certain rights in this software.
 *
 * Redistribution and use in source and binary forms, with or without
 * modification, are permitted provided that the following conditions
 * are met:
 *
 * 1. Redistributions of source code must retain the above copyright
 * notice, this list of conditions and the following disclaimer.
 *
 * 2. Redistributions in binary form must reproduce the above copyright
 * notice, this list of conditions and the following disclaimer in the
 * documentation and/or other materials provided with the distribution.
 *
 * THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
 * "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
 * LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR
 * A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT
 * HOLDER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL,
 * SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT
 * LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE,
 * DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY
 * THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT
 * (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE
 * OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
 */

#include "Predict.h"

#include "BuildFeatures.h"
#include "PredictData.h"

#include <tracktable/RW/KmlOut.h>

using tracktable::kml;

void Predict(const TrajectoryVectorT &_trajectories, const size_t _numSamples) {
    // This routine builds a feature database from specific intervals
    auto features = BuildManyEvenFeatures(_trajectories);
    // use random features as our 'test' set
    auto to_be_predicted = BuildRandomFeatures(_trajectories, 0.2, 0.8);

    // Create rtree for our predictions
    PredictRtreeT rtree;
    // std::vector<PredictRtreeT::value_type> data;

    // Build the feature vector/id number combo for the rtree.  There is an
    // unused value set to 0 that might be used in the future.

    for (auto i = 0u; i < features.size(); ++i) {
        rtree.insert(&(features[i]));
    }

    // Define the the number of trajectories that will be used to predict
    // the destination.  The bins vector will be used to hold the result of
    // how far down the potential list of predicted destinations you need to
    // go to find the right destination.  The bins vector is a little larger
    // to hold the count of getting it wrong.

    std::vector<size_t> bins(_numSamples + 1);

    // Okay.  Here is where the work is done.  Go through each flight and
    // find all of its neighbors to predict where it will land.

    for (auto &current : to_be_predicted) {
        std::vector<PredictRtreeT::value_type> result_n;

        // Note we are getting more results than _numSamples.  This is because
        // we will throw out the hit that corresponds to the trajectory itself.
        // It would be cheating to use that for prediction.

        auto it = rtree.qbegin(boost::geometry::index::nearest(current.feature, _numSamples + 10));
        for (; (it != rtree.qend()) && (result_n.size() < _numSamples); ++it) {
            if ((*it)->index != current.index) result_n.push_back(*it);
        }

        using WeightPairT = std::pair<std::string, double>;
        std::map<WeightPairT::first_type, WeightPairT::second_type> weights;

        TrajectoryVectorT results;

        auto dest = current.index->front().string_property("dest");
        std::cout << dest << std::endl;

        // Take the results from the rtree query, and then build a vector that
        // has the resulting flights.  In addition, build a table of weights for
        // each potential destination (via a map) using what is essentially a
        // 1/d^2 weight.  The d^2 term comes from the "comparable_distance"
        // function.

        auto total_weight = 0.0;
        for (auto &found : result_n) {
            double weight =
                1.0 / (0.01 + boost::geometry::comparable_distance(current.feature, found->feature));
            total_weight += weight;
            results.push_back(found->index);
            weights[found->index->front().string_property("dest")] += weight;
        }

        // An intermediate step.  Basically, the elements of the map are sorted
        // by key (destination) not by value (weight).  We put them in a vector
        // that will be sorted by value since that is how we will use them.

        std::vector<WeightPairT> sorted(weights.begin(), weights.end());

        // Do the actual sorting.  Have to specify using the second element since
        // sort will use the first element (destination) by default.

        std::sort(sorted.begin(), sorted.end(),
                  [](WeightPairT &_lhs, WeightPairT &_rhs) { return _lhs.second > _rhs.second; });

        // Here is where we fill the bin of whether they were ther first guess,
        // second guess, etc., or wrong. The wrong answers go in the 0 bin.
        // The number is done by iterator subtraction, which is a totally
        // valid way to do things in C++.  Note that if pos == sorted.size(),
        // then the find_if returned sorted.end() and the value was
        // not found.

        // size_t cast safe because gauranteed >= 0
        auto pos = size_t(std::find_if(sorted.begin(), sorted.end(),
                                       [&dest](WeightPairT &_wp) { return _wp.first == dest; }) -
                          sorted.begin());
        ++bins[pos == sorted.size() ? 0 : pos + 1];

        if ((pos != sorted.size()) && (pos == 3) && (total_weight > 400)) {
            std::string filename = "output/" + current.index->object_id() + "cand.kml";
            std::ofstream canidateOut(filename);
            if (!canidateOut.is_open()) {
                std::cerr << "Could not open" << filename << std::endl;
            } else {
                canidateOut << kml(results);
            }
            filename = "output/" + current.index->object_id() + ".kml";
            std::ofstream resultOut(filename);
            if (!resultOut.is_open()) {
                std::cerr << "Could not open" << filename << std::endl;
            } else {
                resultOut << kml(*(current.index));
            }
        }
    }

    int total = 0;
    for (unsigned int i = 1; i < bins.size(); ++i) {
        total += bins[i];
        std::cout << "bins[" << i << "] = " << bins[i] << ", total = " << total
                  << ", cumulative fraction  = " << double(total) / double(to_be_predicted.size())
                  << std::endl;
    }
    std::cout << "Got " << bins[0] << " (" << double(bins[0]) / double(to_be_predicted.size())
              << " fraction) wrong" << std::endl;
}

---

PredictData.h

/*
 * Copyright (c) 2014-2023 National Technology and Engineering
 * Solutions of Sandia, LLC. Under the terms of Contract DE-NA0003525
 * with National Technology and Engineering Solutions of Sandia, LLC,
 * the U.S. Government retains certain rights in this software.
 *
 * Redistribution and use in source and binary forms, with or without
 * modification, are permitted provided that the following conditions
 * are met:
 *
 * 1. Redistributions of source code must retain the above copyright
 * notice, this list of conditions and the following disclaimer.
 *
 * 2. Redistributions in binary form must reproduce the above copyright
 * notice, this list of conditions and the following disclaimer in the
 * documentation and/or other materials provided with the distribution.
 *
 * THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
 * "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
 * LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR
 * A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT
 * HOLDER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL,
 * SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT
 * LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE,
 * DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY
 * THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT
 * (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE
 * OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
 */

#ifndef PredictData_h
#define PredictData_h

#include <tracktable/Analysis/GuardedBoostGeometryRTreeHeader.h>
#include <tracktable/Domain/Terrestrial.h>

#include <boost/geometry/geometries/adapted/std_array.hpp>

#include <array>
#include <memory>

using TrajectoryT = tracktable::domain::terrestrial::trajectory_type;
// TODO: Make TrajectoryT a template param

class PredictData {
   public:
    // Have to use boost array instead of standard because rtree can't handle std::array
    using FeatureT = std::array<double, 10>;
    FeatureT feature;
    size_t Id;
    std::shared_ptr<TrajectoryT> index;

    PredictData() = delete;
    PredictData(FeatureT const& _f, size_t id, std::shared_ptr<TrajectoryT> idx)
        : feature(_f), Id(id), index(idx) {}
    ~PredictData() = default;

    PredictData(PredictData const& _other) {
        feature = _other.feature;
        Id = _other.Id;
        index = _other.index;
    }

    PredictData& operator=(PredictData const& _other) {
        feature = _other.feature;
        Id = _other.Id;
        index = _other.index;
        return *this;
    }
};

BOOST_GEOMETRY_REGISTER_STD_ARRAY_CS(cs::cartesian)

namespace boost {
namespace geometry {
namespace index {

template <>
struct indexable<PredictData*> {
    using result_type = PredictData::FeatureT const&;  // required by boost
    result_type& operator()(PredictData* const& v) const { return v->feature; }
};

}  // namespace index
}  // namespace geometry
}  // namespace boost
using PredictRtreeT = boost::geometry::index::rtree<PredictData*, boost::geometry::index::quadratic<16>>;

#endif

---

BuildFeatures.h

/*
 * Copyright (c) 2014-2023 National Technology and Engineering
 * Solutions of Sandia, LLC. Under the terms of Contract DE-NA0003525
 * with National Technology and Engineering Solutions of Sandia, LLC,
 * the U.S. Government retains certain rights in this software.
 *
 * Redistribution and use in source and binary forms, with or without
 * modification, are permitted provided that the following conditions
 * are met:
 *
 * 1. Redistributions of source code must retain the above copyright
 * notice, this list of conditions and the following disclaimer.
 *
 * 2. Redistributions in binary form must reproduce the above copyright
 * notice, this list of conditions and the following disclaimer in the
 * documentation and/or other materials provided with the distribution.
 *
 * THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
 * "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
 * LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR
 * A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT
 * HOLDER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL,
 * SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT
 * LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE,
 * DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY
 * THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT
 * (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE
 * OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
 */
#ifndef BuildFeatures_h
#define BuildFeatures_h

#include "PredictData.h"

#include <tracktable/Domain/Terrestrial.h>

#include <vector>

using TrajectoryT = tracktable::domain::terrestrial::trajectory_type;
using TrajectoryVectorT = std::vector<std::shared_ptr<TrajectoryT>>;
using DataVectorT = std::vector<PredictData>;

DataVectorT BuildFeatures(const TrajectoryVectorT &_trajectories, double _fraction);
DataVectorT BuildManyEvenFeatures(const TrajectoryVectorT &_trajectories);
DataVectorT BuildManyRandomFeatures(const TrajectoryVectorT &_trajectories);
DataVectorT BuildRandomFeatures(const TrajectoryVectorT &_trajectories, double _lower, double _upper);
PredictData BuildFeature(std::shared_ptr<TrajectoryT> _trajectory, double _fraction);

#endif

---

BuildFeatures.cpp

/*
 * Copyright (c) 2014-2023 National Technology and Engineering
 * Solutions of Sandia, LLC. Under the terms of Contract DE-NA0003525
 * with National Technology and Engineering Solutions of Sandia, LLC,
 * the U.S. Government retains certain rights in this software.
 *
 * Redistribution and use in source and binary forms, with or without
 * modification, are permitted provided that the following conditions
 * are met:
 *
 * 1. Redistributions of source code must retain the above copyright
 * notice, this list of conditions and the following disclaimer.
 *
 * 2. Redistributions in binary form must reproduce the above copyright
 * notice, this list of conditions and the following disclaimer in the
 * documentation and/or other materials provided with the distribution.
 *
 * THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
 * "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
 * LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR
 * A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT
 * HOLDER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL,
 * SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT
 * LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE,
 * DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY
 * THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT
 * (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE
 * OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
 */

#include "BuildFeatures.h"

#include <algorithm>  // for transform
#include <iterator>   //for back_insterter
#include <random>

DataVectorT BuildFeatures(const TrajectoryVectorT &_trajectories, double _fraction) {
    DataVectorT features;
    // transform allows creating a vector of type b from a vector of type a
    std::transform(_trajectories.begin(), _trajectories.end(), std::back_inserter(features),
                   [&](std::shared_ptr<TrajectoryT> _t) { return BuildFeature(_t, _fraction); });
    return features;
}

DataVectorT BuildManyEvenFeatures(const TrajectoryVectorT &_trajectories) {
    DataVectorT features;
    for (auto i = 2u; i <= 8u; ++i) {
        auto f = BuildFeatures(_trajectories, i / 10.0);
        features.insert(features.end(), f.begin(), f.end());
    }
    return features;
}

DataVectorT BuildManyRandomFeatures(const TrajectoryVectorT &_trajectories) {
    DataVectorT features;
    for (auto i = 0u; i < 7u; ++i) {
        auto f = BuildRandomFeatures(_trajectories, 0.2, 0.8);
        features.insert(features.end(), f.begin(), f.end());
    }
    return features;
}

DataVectorT BuildRandomFeatures(const TrajectoryVectorT &_trajectories, double _lower, double _upper) {
    static std::mt19937 gen(0);
    DataVectorT features;
    std::uniform_real_distribution<double> dis(_lower, _upper);  // uniform distribution
    std::transform(_trajectories.begin(), _trajectories.end(), std::back_inserter(features),
                   [&](std::shared_ptr<TrajectoryT> _t) { return BuildFeature(_t, dis(gen)); });
    return features;
}

PredictData BuildFeature(std::shared_ptr<TrajectoryT> _trajectory, double _fraction) {
    // Nothing complicated here.  Basically, divide the trajectory up into
    // some even-sized pieces, and take the early way-points and build a
    // feature vector out of them. The only tricky part is that when you
    // divide a trajectory into n pieces, you essentially have n+1 points
    // to choose from.

    constexpr auto numSamples = 4u;  // Must be >= 2 and <= 4

    PredictData::FeatureT feature;
    static_assert(((numSamples - 1) * 2 + 1) < feature.size(), "Size Mismatch");
    for (auto i = 0u; i < numSamples; ++i) {
        auto frac = _fraction * double(i) / double(numSamples - 1.0);
        auto p = tracktable::point_at_length_fraction(*_trajectory, frac);
        feature[2 * i] = p.longitude();
        feature[2 * i + 1] = p.latitude();
    }
    auto p = tracktable::point_at_length_fraction(*_trajectory, _fraction);
    feature[2 * numSamples] = (p.timestamp() - _trajectory->front().timestamp()).seconds() / 1200.;

    return PredictData(feature, 0, _trajectory);
}