landscape.hpp

// Copyright 2021 Arnaud Becheler    <abechele@umich.edu>
 
/*********************************************************************** * This program is free software; you can
 *redistribute it and/or modify * it under the terms of the GNU General Public License as published by * the Free
 *Software Foundation; either version 2 of the License, or    * (at your option) any later version. *
 *                                                                      *
 ***************************************************************************/
 
#pragma once
 
#include "raster.hpp"
 
#include <algorithm>
#include <assert.h>
#include <functional> // std::cref
#include <iostream>
#include <map>
#include <optional>
#include <set>
#include <stdexcept>
#include <string>
#include <tuple> // std::forward_as_tuple
 
namespace quetzal::geography
{
template <typename Key = std::string, typename Time = int> class landscape
{
 
  public:
    using key_type = Key;
 
    using time_type = Time;
 
    using location_descriptor = typename raster<time_type>::location_descriptor;
 
    using time_descriptor = typename raster<time_type>::time_descriptor;
 
    using latlon = typename raster<time_type>::latlon;
 
    using lonlat = typename raster<time_type>::lonlat;
 
    using colrow = typename raster<time_type>::colrow;
 
    using rowcol = typename raster<time_type>::rowcol;
 
    using decimal_degree = typename raster<time_type>::decimal_degree;
 
  private:
    std::map<key_type, raster<Time>> _variables;
 
    landscape(const std::map<key_type, std::filesystem::path> &files, const std::vector<time_type> &times)
        : _variables()
    {
 
        for (auto const &it : files)
        {
            _variables.emplace(std::piecewise_construct, std::forward_as_tuple(it.first),
                               std::forward_as_tuple(it.second, times));
        }
 
        auto have_same_origin = [this](auto const &it) {
            return it.second.origin() == _variables.cbegin()->second.origin();
        };
 
        if (!all_of(_variables.cbegin(), _variables.cend(), have_same_origin))
        {
            throw std::runtime_error("all ecological quantities dataset must have same origin");
        }
 
        auto have_same_extent = [this](auto const &it) {
            return it.second.get_extent() == _variables.cbegin()->second.get_extent();
        };
 
        if (!all_of(_variables.cbegin(), _variables.cend(), have_same_extent))
        {
            throw std::runtime_error("all ecological quantities dataset must have same extent");
        }
 
        auto have_same_resolution = [this](auto const &it) {
            return it.second.get_resolution() == _variables.cbegin()->second.get_resolution();
        };
 
        if (!all_of(_variables.cbegin(), _variables.cend(), have_same_resolution))
        {
            throw std::runtime_error("all ecological quantities dataset must have same resolution");
        }
 
        auto have_same_depth = [this](auto const &it) {
            return it.second.depth() == _variables.cbegin()->second.depth();
        };
 
        if (!all_of(_variables.cbegin(), _variables.cend(), have_same_depth))
        {
            throw std::runtime_error("all ecological quantities dataset must have same depth");
        }
    }
 
  public:
 
    inline static landscape from_file(const std::map<key_type, std::filesystem::path> &files,
                                const std::vector<time_type> &times)
    {
        return landscape(files, times);
    }
 
    template <typename Callable>
        requires std::invocable<Callable, location_descriptor, time_descriptor, std::optional<double>>
    inline void to_geotiff(Callable f, time_descriptor start, time_descriptor end, const std::filesystem::path &file) const
    {
        _variables.cbegin()->second.to_geotiff(f, start, end, file);
    }
 
    inline void to_shapefile(std::vector<latlon> points, const std::filesystem::path &file) const
    {
        _variables.cbegin()->second.to_shapefile(points, file);
    }
 
    inline void to_shapefile(std::map<latlon, int> counts, const std::filesystem::path &file) const
    {
        _variables.cbegin()->second.to_shapefile(counts, file);
    }
 
    std::ostream &write(std::ostream &stream) const
    {
        stream << "Landscape of " << _variables.size() << " aligned rasters:\n";
 
        for (const auto &it : _variables)
            stream << it.first << " ";
 
        stream << "\nOrigin: " << origin() << "\nWidth: " << width() << "\nHeight: " << height()
               << "\nDepth: " << depth() << "\nResolution: "
               << "\n\tLat: " << get_resolution().lat() << "\n\tLon: " << get_resolution().lon() << "\nExtent:"
               << "\n\tLat min: " << get_extent().lat_min() << "\n\tLat max: " << get_extent().lat_max()
               << "\n\tLon min: " << get_extent().lon_min() << "\n\tLon max: " << get_extent().lon_max();
        return stream;
    }
 
 
 
    inline auto num_variables() const noexcept
    {
        return _variables.size();
    }
 
    inline auto num_locations() const noexcept
    {
        return this->width() * this->height();
    }
 
 
 
    inline const raster<time_type> &operator[](const key_type &key) const
    {
        return _variables.at(key);
    }
 
    inline raster<time_type> &operator[](const key_type &key)
    {
        return _variables.at(key);
    }
 
 
 
    inline latlon origin() const noexcept
    {
        return _variables.cbegin()->second.origin();
    }
 
    inline resolution<decimal_degree> get_resolution() const noexcept
    {
        return _variables.cbegin()->second.get_resolution();
    }
 
    inline extent<decimal_degree> get_extent() const noexcept
    {
        return _variables.cbegin()->second.get_extent();
    }
 
    inline int width() const noexcept
    {
        return _variables.cbegin()->second.width();
    }
 
    inline int height() const noexcept
    {
        return _variables.cbegin()->second.height();
    }
 
    inline int depth() const noexcept
    {
        return _variables.cbegin()->second.depth();
    }
 
 
 
    inline auto locations() const noexcept
    {
        return _variables.cbegin()->second.locations();
    }
 
    inline auto times() const noexcept
    {
        return _variables.cbegin()->second.times();
    }
 
    inline bool is_valid(location_descriptor x) const noexcept
    {
        return _variables.cbegin()->second.is_valid(x);
    }
 
    inline bool is_valid(const colrow &x) const noexcept
    {
        return _variables.cbegin()->second.is_valid(x);
    }
 
    inline bool is_valid(const rowcol &x) const noexcept
    {
        return _variables.cbegin()->second.is_valid(x);
    }
 
    inline bool contains(const latlon &x) const noexcept
    {
        return _variables.cbegin()->second.contains(x);
    }
 
    inline bool contains(const lonlat &x) const noexcept
    {
        return _variables.cbegin()->second.contains(x);
    }
 
    inline bool is_valid(time_descriptor t) const noexcept
    {
        return _variables.cbegin()->second.is_valid(t);
    }
 
    inline bool is_recorded(const time_type &t) const noexcept
    {
        return _variables.cbegin()->second.is_recorded(t);
    }
 
    inline bool is_in_interval(const time_type &t) const noexcept
    {
        return _variables.cbegin()->second.is_in_interval(t);
    }
 
 
 
    inline location_descriptor to_descriptor(const colrow &x) const noexcept
    {
        assert( is_valid(x) );
        return _variables.cbegin()->second.to_descriptor(x);
    }
 
    inline location_descriptor to_descriptor(const latlon &x) const noexcept
    {
        assert( contains(x) );
        return _variables.cbegin()->second.to_descriptor(x);
    }
 
    colrow to_colrow(const latlon &x) const noexcept
    {
        assert( contains(x) );
        return _variables.cbegin()->second.to_colrow(x);
    }
 
    inline rowcol to_rowcol(const latlon &x) const noexcept
    {
        assert( contains(x) );
        return _variables.cbegin()->second.to_rowcol(x);
    }
 
    inline colrow to_colrow(location_descriptor x) const noexcept
    {
        assert( is_valid(x) );
        return _variables.cbegin()->second.to_colrow(x);
    }
 
    inline latlon to_latlon(location_descriptor x) const noexcept
    {
        assert( is_valid(x) );
        return _variables.cbegin()->second.to_latlon(x);
    }
 
    inline lonlat to_lonlat(location_descriptor x) const noexcept
    {
        assert( is_valid(x) );
        return _variables.cbegin()->second.to_lonlat(x);
    }
    
    inline latlon to_latlon(const colrow &x) const noexcept
    {
        assert( is_valid(x) );
        return _variables.cbegin()->second.to_latlon(x);
    }
 
    inline lonlat to_lonlat(const colrow &x) const noexcept
    {
        assert( is_valid(x) );
        return _variables.cbegin()->second.to_lonlat(x);
    }
 
    inline latlon to_centroid(const latlon &x) const
    {
        assert(this->contains(x));
        return _variables.cbegin()->second.to_centroid(x);
    }
 
};
 
template <typename Key, typename Time> std::ostream &operator<<(std::ostream &os, const landscape<Key, Time> &l)
{
    return l.write(os);
}
 
} // namespace quetzal::geography