Merge pull request #28 from Quetzal-framework/develop

feat: dispersal kernel refactoring, examples and documentation

.vscode/tasks.json

@@ -134,5 +134,17 @@
             "dependsOn": ["Build Template Project (Debug)"],
             "problemMatcher": [],
         },
+        // DOCUMENTATION
+        {
+            "label": "Make Docs",
+            "type": "shell",
+            "command": "cd ${workspaceFolder}/build && make docs ",
+            "group": "test",
+            "presentation": {
+                "reveal": "always",
+                "panel": "shared"
+            },
+            "dependsOn": ["Run Unit Tests (Debug)"],
+        },
     ]
 }

conanfile.py

@@ -42,10 +42,12 @@ def requirements(self):
                 self.requires("boost/[>1.75 <1.80]")
                 self.requires("gdal/[>=3.4.0]")
                 self.requires("range-v3/0.12.0")
+                self.requires("mp-units/2.0.0")
        if self.settings.os == "Linux":
             if self.settings.compiler == "gcc":
                 self.requires("boost/[>1.75 <1.80]")
                 self.requires("gdal/[>=3.4.0 <3.5.1]")
+                self.requires("mp-units/2.0.0")
                 self.requires("range-v3/0.12.0")
                 self.requires("libtiff/4.5.1", override=True) # Version conflict: libgeotiff/1.7.1->libtiff/4.6.0, gdal/3.4.3->libtiff/4.5.1.
                 self.requires("libdeflate/1.18", override=True) # Version conflict: libtiff/4.6.0->libdeflate/1.19, gdal/3.4.3->libdeflate/1.18.

docs/4-tutorials.md

@@ -811,3 +811,30 @@ There a different modeling approaches to look at these demographic histories, pr
 - **Spatial graph structure:** Another category of models puts an emphasis on the geographic structure of these populations and embed GIS information into the model. Here the demographic history receives the semantic of a dense spatial graph where geolocalized nodes represent demes and local processes whereas edges represent distances and dispersal processes. In this kind of models estimating the individual properties of so many demes is not judicious and genetic inference rather aims at estimating species-wide parameters that link local deme conditions (such as rainfall) to populational processes (like growth rate). It's crucial to understand that within this geospatial framework, events like population changes and admixture events are emergent properties of the model: they are generally not explicitly outlined in the model nor subjected to direct estimation.
 
 ---
+
+[//]: # (----------------------------------------------------------------------)
+@page dispersal_kernels Dispersal Kernels
+
+@note 
+The objective of this section is to parametrize a Dispersal Location Kernel (in the sense of <a href="pdf/nathan-et-al-2012.pdf" target="_blank"><b> Nathan et al. 2012</b></a>)
+and to compute useful quantities, such as the distance between two locations, the probability to disperse from one to the other, and the mean dispersal distance expected under the distribution.
+
+## Background 
+
+The dispersal location kernel represents the statistical pattern of dispersal distances within a population. In this context, it essentially serves as the probability density function (pdf) that outlines the distribution of locations after dispersal in relation to the original location. The expressions have been adjusted to incorporate a scale parameter, 
+\f$a\f$, which is consistent with a distance unit, and a shape parameter, 
+\f$b\f$, that dictates the curvature of the distribution curve, emphasizing the influence of long-distance events.
+
+For a source \f$(x_0,y_0)\f$, the dispersal location kernel denoted as \f$k_L(r)\f$ provides the density of the probability of the dispersal end point in the 2D space. In this case, \f$k_L(r)dA\f$ is the probability of a dispersal end point to be within a small 2D area \f$dA\f$ around the location \f$(x,y)\f$. Since a probability is unitless and \f$dA\f$ is an area, \f$k_L(r)\f$ is expressed in per unit area in a 2D space.
+
+Quetzal automate dimensional analysis and conversion thanks to the `mp-units` library.
+
+**Input**
+
+@include{lineno} dispersal_kernel.cpp
+
+**Output**
+
+@include{lineno} dispersal_kernel.txt
+
+---

docs/CMakeLists.txt

@@ -17,7 +17,7 @@ else ()
   # doxygen settings can be set here, prefixed with "DOXYGEN_"
   set(DOXYGEN_PROJECT_NAME "Quetzal-CoaTL")
   set(DOXYGEN_PROJECT_BRIEF "The Coalescence Template Library")
-  set(DOXYGEN_PROJECT_LOGO "quetzal-coaltl-fa.png")
+  set(DOXYGEN_PROJECT_LOGO "assets/images/quetzal-coaltl-fa.png")
   set(DOXYGEN_USE_MDFILE_AS_MAINPAGE "mainpage.md")
   set(DOXYGEN_EXCLUDE_PATTERNS "README.md" "doxygen-awesome-css" "*env*" "conanfile.*")
   set(DOXYGEN_OUTPUT_DIRECTORY "${PROJECT_BINARY_DIR}/docs")
@@ -28,7 +28,8 @@ else ()
                                     "doxygen-awesome-css/doxygen-awesome-fragment-copy-button.js")
   # Documentation for tutorials requires ctest call to populate build folder
   set(DOXYGEN_EXAMPLE_PATH "../example" "../build/example/output")
-  set(DOXYGEN_IMAGE_PATH "assets")
+  set(DOXYGEN_IMAGE_PATH "assets/images")
+  set(DOXYGEN_HTML_EXTRA_FILES = "assets/pdf")
 
   # Required by doxygen awesome
   set(DOXYGEN_GENERATE_TREEVIEW YES)

docs/mainpage.md

@@ -22,4 +22,4 @@ encourage reuse of algorithms and data structures.
 
 This library intends to do this for the C++ language.
 
-\image html pictures/quetzal_scheme.png
+\image html assets/images/quetzal-scheme.png

example/CMakeLists.txt

@@ -10,6 +10,9 @@ find_package(GDAL REQUIRED)
 # range-v3
 find_package(range-v3)
 
+# mp-units
+find_package(mp-units CONFIG REQUIRED)
+
 # Include Boost as an imported target: 1.73 important for compilation with gcc 10 and C++20
 find_package(Boost 1.79 REQUIRED COMPONENTS unit_test_framework filesystem serialization)
 if(Boost_FOUND)
@@ -33,7 +36,7 @@ foreach(testSrc ${SRCS})
   # Ignore warnings about subtle ABI change
   target_compile_options(${testName} PUBLIC "-Wno-psabi")
   # Link to Boost libraries AND other targets and dependencies
-  target_link_libraries(${testName} quetzal::quetzal boost::boost GDAL::GDAL range-v3::range-v3)
+  target_link_libraries(${testName} quetzal::quetzal boost::boost GDAL::GDAL range-v3::range-v3 mp-units::mp-units)
   # Specifies include directories to use when compiling a given target
   target_include_directories(
     ${testName} PRIVATE

example/dispersal_kernel.cpp

@@ -0,0 +1,30 @@
+#include "quetzal/quetzal.hpp"
+#include <iostream>
+#include <mp-units/ostream.h>
+
+using namespace quetzal;
+
+int main()
+{
+	// Shorten notation
+	using namespace mp_units::si::unit_symbols; // SI system: km, m, s 
+	using quetzal::demography::dispersal_kernel::Gaussian;
+	using quetzal::geography::lonlat;
+
+	// Dispersal kernels operate on distances between geographic coordinates
+    constexpr auto source = lonlat(2.25, 48.9176);
+    constexpr auto target = lonlat(2.25, 48.95);
+	constexpr auto distance = source.great_circle_distance_to(target) * km;
+
+	// Used to parametrize the kernel
+	constexpr auto param = Gaussian<>::param_type {.a=1000.*m};
+
+	// Compute quantities
+	auto p = Gaussian<>::pdf(distance, param);
+	auto d = Gaussian<>::mean_dispersal_distance(param);
+
+	std::cout << "Dispersal from " << source << " to " << target << " :\n"
+			  << "distance is " << distance << "\n"
+			  << "probability is " << p << "\n"
+			  << "mean dispersal distance is " << d << std::endl; 
+}

example/geography_raster_1.cpp

@@ -30,7 +30,7 @@ using namespace quetzal;
 	//
 	//          90
 	//          |         (+)
-	// 	    0 --------------> 180    X size positive in decimal degre (east direction positive)
+	// 	    0 --------------> 180    X size positive in decimal degree (east direction positive)
 	//          |                    Y size negative in decimal degree (south direction negative)
 	//          |
 	// 			|

example/spatial_graph.cpp

@@ -0,0 +1,32 @@
+#include "quetzal/quetzal.hpp"
+#include <filesystem>
+#include <cassert>
+
+using namespace quetzal;
+
+int main()
+{
+	auto file1 = std::filesystem::current_path() / "data/bio1.tif";
+	auto file2 = std::filesystem::current_path() / "data/bio12.tif";
+
+	// The raster have 10 bands that we will assign to 2001 ... 2010.
+	std::vector<int> times(10);
+    std::iota(times.begin(), times.end(), 2001);
+
+	// Initialize the landscape: for each var a key and a file, for all a time series.
+	using landscape_type = quetzal::geography::landscape<>;
+	auto env = quetzal::geography::landscape<>::from_files( { {"bio1", file1}, {"bio12", file2} }, times );
+	std::cout << env << std::endl;
+
+	// // Embed default demographic quantities along vertices (N_{xt}) and edges (Phi_{xyt})
+	// auto graph1 = env.to_4_neighborhood_graph<>();
+	// auto graph2 = env.to_8_neighborhood_graph<>();
+	// auto graph3 = env.to_complete_graph<>();
+
+	// auto weights graph.edges()
+	// 	| [](const & e){ return { graph.source(e), graph.target(e)};  }
+	// 	| [](const & p){ return { env.to_latlon(s), graph.to_latlon(t)};  }
+	// 	| [](const & p){ return x.great_circle_distance_to(y);}
+	// 	| [](const & d){ return quetzal::demography::dispersal_kernel::ExponentialPower(d, {.a=1., .b=5.5});}
+	// 	| []()
+}