Main Lemma Repository
You can not select more than 25 topics Topics must start with a letter or number, can include dashes ('-') and can be up to 35 characters long.

pyLemmaCore.cpp 7.2KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134
  1. /* This file is part of Lemma, a geophysical modelling and inversion API.
  2. * More information is available at http://lemmasoftware.org
  3. */
  4. /* This Source Code Form is subject to the terms of the Mozilla Public
  5. * License, v. 2.0. If a copy of the MPL was not distributed with this
  6. * file, You can obtain one at http://mozilla.org/MPL/2.0/.
  7. */
  8. /**
  9. * @file
  10. * @date 22/04/19 14:06:32
  11. * @version $Id$
  12. * @author Trevor Irons (ti)
  13. * @email Trevor.Irons@utah.edu
  14. * @copyright Copyright (c) 2019, University of Utah
  15. * @copyright Copyright (c) 2019, Lemma Software, LLC
  16. */
  17. #include <pybind11/pybind11.h>
  18. #include <pybind11/iostream.h>
  19. #include <pybind11/eigen.h>
  20. #include "LemmaCore"
  21. namespace py = pybind11;
  22. PYBIND11_MODULE(LemmaCore, m) {
  23. py::add_ostream_redirect(m, "ostream_redirect");
  24. m.doc() = "Python binding of LemmaCore, additional details can be found at https://lemmasoftware.org";
  25. py::class_<Lemma::RectilinearGrid, std::shared_ptr<Lemma::RectilinearGrid> > RectilinearGrid(m, "RectilinearGrid");
  26. // lifecycle
  27. RectilinearGrid.def(py::init(&Lemma::RectilinearGrid::NewSP))
  28. .def_static("DeSerialize", py::overload_cast<const std::string&>(&Lemma::RectilinearGrid::DeSerialize), "Construct object from yaml representation")
  29. // print
  30. .def("__repr__", &Lemma::RectilinearGrid::Print)
  31. .def("Serialize", &Lemma::RectilinearGrid::Print, "YAML representation of the class")
  32. // accessors
  33. .def("GetName", &Lemma::RectilinearGrid::GetName, "Returns the name of the class")
  34. .def("GetNx", &Lemma::RectilinearGrid::GetNx, "Returns the number of cells in the x direction")
  35. .def("GetNy", &Lemma::RectilinearGrid::GetNy, "Returns the number of cells in the y direction")
  36. .def("GetNz", &Lemma::RectilinearGrid::GetNz, "Returns the number of cells in the z direction")
  37. .def("GetOx", &Lemma::RectilinearGrid::GetOx, "Returns the grid origin offset in the x direction")
  38. .def("GetOy", &Lemma::RectilinearGrid::GetOy, "Returns the grid origin offset in the y direction")
  39. .def("GetOz", &Lemma::RectilinearGrid::GetOz, "Returns the grid origin offset in the z direction")
  40. .def("GetDx", &Lemma::RectilinearGrid::GetDx, "Returns the grid spacing in the x direction")
  41. .def("GetDy", &Lemma::RectilinearGrid::GetDy, "Returns the grid spacing in the y direction")
  42. .def("GetDz", &Lemma::RectilinearGrid::GetDz, "Returns the grid spacing in the z direction")
  43. // modifiers
  44. .def("SetDimensions", &Lemma::RectilinearGrid::SetDimensions, "Sets the number of cells in x, y, and z")
  45. .def("SetOffset", &Lemma::RectilinearGrid::SetOffset, "Sets the origin offset in x, y, and z")
  46. .def("SetSpacing", &Lemma::RectilinearGrid::SetSpacing, "Sets the grid spacing in x, y, and z")
  47. ;
  48. py::class_<Lemma::ASCIIParser, std::shared_ptr<Lemma::ASCIIParser> > ASCIIParser(m, "ASCIIParser");
  49. // lifecycle
  50. ASCIIParser.def(py::init(&Lemma::ASCIIParser::NewSP))
  51. .def_static("DeSerialize", py::overload_cast<const std::string&>(&Lemma::ASCIIParser::DeSerialize), "Construct object from yaml representation")
  52. // print
  53. .def("__repr__", &Lemma::ASCIIParser::Print)
  54. .def("Serialize", &Lemma::ASCIIParser::Print, "YAML representation of the class")
  55. // accessors
  56. .def("GetName", &Lemma::ASCIIParser::GetName, "Returns the name of the class")
  57. .def("GetFileLocation", &Lemma::ASCIIParser::GetFileLocation, "Returns the current file location")
  58. // modifiers
  59. .def("SetCommentString", &Lemma::ASCIIParser::SetCommentString, "Sets the comment string after which all text is ignored")
  60. .def("SetBufferSize", &Lemma::ASCIIParser::SetBufferSize, "Sets the buffer size")
  61. // methods
  62. .def("Open", &Lemma::ASCIIParser::Open, "Opens file specified by argument")
  63. .def("Close", &Lemma::ASCIIParser::Close, "Closes current file object")
  64. .def("ReadReals", &Lemma::ASCIIParser::ReadReals, "Returns vector of nr reals")
  65. .def("ReadInts", &Lemma::ASCIIParser::ReadInts, "Returns vector of ni ints")
  66. .def("ReadStrings", &Lemma::ASCIIParser::ReadStrings, "Returns vector of ns strings")
  67. .def("JumpToLocation", &Lemma::ASCIIParser::JumpToLocation, "File object jumps to specified location")
  68. ;
  69. py::class_<Lemma::CubicSplineInterpolator, std::shared_ptr<Lemma::CubicSplineInterpolator> >(m, "CubicSplineInterpolator")
  70. // lifecycle
  71. .def(py::init(&Lemma::CubicSplineInterpolator::NewSP))
  72. .def_static("DeSerialize", py::overload_cast<const std::string&>(&Lemma::CubicSplineInterpolator::DeSerialize), "Construct object from yaml representation")
  73. // print
  74. .def("__repr__", &Lemma::CubicSplineInterpolator::Print)
  75. .def("Serialize", &Lemma::CubicSplineInterpolator::Print, "YAML representation of the class")
  76. // accessors
  77. .def("GetName", &Lemma::CubicSplineInterpolator::GetName, "Returns the name of the class")
  78. .def("GetKnotAbscissa", &Lemma::CubicSplineInterpolator::GetKnotAbscissa, "Returns the knot abscissa values")
  79. .def("GetKnotOrdinate", &Lemma::CubicSplineInterpolator::GetKnotOrdinate, "Returns the knot ordinate values")
  80. // modifiers
  81. .def("SetKnots", &Lemma::CubicSplineInterpolator::SetKnots, "Sets the knots to use for interpolation")
  82. .def("ResetKnotOrdinate", &Lemma::CubicSplineInterpolator::ResetKnotOrdinate, "Resets the knots to use for interpolation, when abscissa values haven't changed")
  83. // methods
  84. .def("InterpolateOrderedSet", &Lemma::CubicSplineInterpolator::InterpolateOrderedSet, "Interpolate a monotonically increasing ordered set.")
  85. .def("Integrate", py::overload_cast<const Lemma::Real&, const Lemma::Real& >(&Lemma::CubicSplineInterpolator::Integrate), "Integrates between the arguments using cubic spline values.")
  86. .def("IntegrateN", py::overload_cast<const Lemma::Real&, const Lemma::Real&, const int& >(&Lemma::CubicSplineInterpolator::Integrate), "Integrates the spline from x0 to x1. Uses composite Simpson's rule and n is the number of segments")
  87. .def("Interpolate", py::overload_cast<const Lemma::Real& >(&Lemma::CubicSplineInterpolator::Interpolate), "Interpolation at a single point, x is the interpolation abscissa point, returns the ordinate value at x")
  88. .def("InterpolateI", py::overload_cast<const Lemma::Real&, int& >(&Lemma::CubicSplineInterpolator::Interpolate), "Interpolation at a single point, x is the interpolation abscissa point and i is the knot to begin searchin at returns the ordinate value at x")
  89. ;
  90. // ABC
  91. //py::class_<Lemma::EarthModel, std::shared_ptr<Lemma::EarthModel> >(m, "EarthModel")
  92. /*
  93. py::class_<Lemma::LayeredEarth, std::shared_ptr<Lemma::LayeredEarth> >(m, "LayeredEarth")
  94. // lifecycle
  95. .def(py::init(&Lemma::LayeredEarth::NewSP))
  96. //.def_static("DeSerialize", py::overload_cast<const std::string&>(&Lemma::LayeredEarth::DeSerialize), "Construct object from yaml representation")
  97. // print
  98. .def("__repr__", &Lemma::LayeredEarth::Print)
  99. .def("Serialize", &Lemma::LayeredEarth::Print, "YAML representation of the class")
  100. ;
  101. */
  102. }