lsst.meas.extensions.astrometryNet  master-g43362ee2f3+22
 All Classes Namespaces Files Functions Variables Properties Groups Pages
astrometry_net.cc
Go to the documentation of this file.
1 // -*- lsst-C++ -*-
2 
3 // Astrometry.net include files...
4 extern "C" {
5 #include "astrometry/solver.h"
6 #include "astrometry/index.h"
7 #include "astrometry/multiindex.h"
8 #include "astrometry/starkd.h"
9 #include "astrometry/fitsioutils.h"
10 #include "astrometry/fitstable.h"
11 #include "astrometry/log.h"
12 #include "astrometry/tic.h"
13 #include "astrometry/healpix.h"
14 
15 #undef ATTRIB_FORMAT
16 #undef FALSE
17 #undef TRUE
18 
19 #undef logdebug
20 #undef debug
21 }
22 
23 #include <string>
24 #include <sstream>
25 
26 #include "pybind11/pybind11.h"
27 #include "pybind11/stl.h"
28 
29 #include "lsst/log/Log.h"
30 #include "lsst/utils/python.h"
32 
33 namespace py = pybind11;
34 using namespace pybind11::literals;
35 
36 namespace lsst {
37 namespace meas {
38 namespace extensions {
39 namespace astrometryNet {
40 namespace {
41 
42 /*
43  * Wrap index_t
44  *
45  * This class has no constructor; the wrapper is only used for objects returned by MultiIndex::operator[]
46  */
47 static void declareIndex(py::module& mod) {
48  py::class_<index_t> cls(mod, "index_t");
49 
50  cls.def("overlapsScaleRange",
51  [](index_t& self, double qlo, double qhi) { return index_overlaps_scale_range(&self, qlo, qhi); },
52  "qlow"_a, "qhigh"_a);
53  cls.def("reload", [](index_t& self) {
54  if (index_reload(&self)) {
55  std::ostringstream os;
56  os << "Failed to reload multi-index file " << self.indexname;
57  throw LSST_EXCEPT(lsst::pex::exceptions::RuntimeError, os.str());
58  }
59  });
60 
61  cls.def_readonly("indexid", &index_t::indexid);
62  cls.def_readonly("healpix", &index_t::healpix);
63  cls.def_readonly("hpnside", &index_t::hpnside);
64  cls.def_readonly("nstars", &index_t::nstars);
65  cls.def_readonly("nquads", &index_t::nquads);
66 }
67 
71 static void declareMultiIndex(py::module& mod) {
72  py::class_<MultiIndex> cls(mod, "MultiIndex");
73 
74  cls.def(py::init<std::string const&>(), "filepath"_a);
75 
76  cls.def("__getitem__",
77  [](MultiIndex const& self, int i) {
78  auto cind = lsst::utils::python::cppIndex(self.getLength(), i);
79  return self[cind];
80  },
81  py::return_value_policy::reference_internal, py::is_operator());
82 
83  cls.def("addIndex", &MultiIndex::addIndex, "filepath"_a, "metadataOnly"_a);
84  cls.def("isWithinRange", &MultiIndex::isWithinRange, "ra"_a, "dec"_a, "radius"_a);
85  cls.def("unload", &MultiIndex::unload);
86  cls.def_property_readonly("name", &MultiIndex::getName);
87  cls.def("__len__", &MultiIndex::getLength);
88  cls.def("reload", &MultiIndex::reload);
89 }
90 
94 static void declareSolver(py::module& mod) {
95  py::class_<Solver> cls(mod, "Solver");
96 
97  cls.def(py::init<>());
98 
126  cls.def("getCatalog", &Solver::getCatalog, "inds"_a, "ctrCoord"_a, "radius"_a, "idCol"_a,
127  "filterNameList"_a, "magColList"_a, "magErrColList"_a, "starGalCol"_a, "varCol"_a,
128  "uniqueIds"_a = true);
129  cls.def("getSolveStats", &Solver::getSolveStats);
130  cls.def("getWcs", &Solver::getWcs);
131  cls.def("didSolve", &Solver::didSolve);
132  cls.def("run", &Solver::run, "cpulimit"_a);
133  cls.def("getQuadSizeRangeArcsec", &Solver::getQuadSizeRangeArcsec);
134  cls.def("addIndices", &Solver::addIndices, "indices"_a);
135  cls.def("setParity", &Solver::setParity, "setParityFlipped", "parity"_a);
136  cls.def("setMatchThreshold", &Solver::setMatchThreshold, "threshold"_a);
137  cls.def("setPixelScaleRange", &Solver::setPixelScaleRange, "low"_a, "high"_a);
138  cls.def("setRaDecRadius", &Solver::setRaDecRadius, "ra"_a, "dec"_a, "rad"_a);
139  cls.def("setImageSize", &Solver::setImageSize, "width"_a, "height"_a);
140  cls.def("setMaxStars", &Solver::setMaxStars, "maxStars"_a);
141  cls.def("setStars", &Solver::setStars, "sourceCat"_a, "x0"_a, "y0"_a);
142 }
143 
144 // declare logging functions for use by the Python
145 
146 LOG_LOGGER an_log = LOG_GET("meas.astrom.astrometry_net");
147 
148 void an_log_callback(void* baton, enum log_level level, const char* file, int line, const char* func,
149  const char* format, va_list va) {
150  // translate between logging levels
151  int levelmap[5];
152  levelmap[LOG_NONE] = LOG_LVL_FATAL;
153  levelmap[LOG_ERROR] = LOG_LVL_FATAL;
154  levelmap[LOG_MSG] = LOG_LVL_INFO;
155  levelmap[LOG_VERB] = LOG_LVL_DEBUG;
156  levelmap[LOG_ALL] = LOG_LVL_DEBUG;
157  int lsstlevel = levelmap[level];
158 
159  va_list vb;
160  // find out how long the formatted string will be
161  va_copy(vb, va);
162  const int len = vsnprintf(NULL, 0, format, va) + 1; // "+ 1" for the '\0'
163  va_end(va);
164  // allocate a string of the appropriate length
165  char msg[len];
166  (void)vsnprintf(msg, len, format, vb);
167  va_end(vb);
168 
169  // trim trailing \n
170  if (msg[len - 2] == '\n') {
171  msg[len - 2] = '\0';
172  }
173 
174  // not using the LOGS macro because the original location info is wanted
175  if (an_log.isEnabledFor(lsstlevel)) {
176  an_log.logMsg(log4cxx::Level::toLevel(lsstlevel), log4cxx::spi::LocationInfo(file, func, line), msg);
177  }
178 }
179 
181 void start_an_logging() {
182  // NOTE, this has to happen before the log_use_function!
183  log_init(LOG_VERB);
184  log_use_function(an_log_callback, NULL);
185  log_to(NULL);
186 }
187 
189 void finalize() {
190  log_use_function(NULL, NULL);
191  log_to(stdout);
192 }
193 
194 } // namespace <anonymous>
195 
196 PYBIND11_PLUGIN(astrometry_net) {
197  py::module mod("astrometry_net");
198 
199  // code that is run at import time
200  fits_use_error_system();
201  start_an_logging();
202 
203  mod.def("healpixDistance", &healpixDistance, "hp"_a, "nside"_a, "coord"_a);
204 
205  mod.def("an_log_init", [](int level) { log_init(static_cast<log_level>(level)); }, "level"_a);
206 
207  mod.def("an_log_get_level", []() { return static_cast<int>(log_get_level()); });
208  mod.def("an_log_set_level", [](int level) { log_set_level(static_cast<log_level>(level)); }, "level"_a);
209  mod.def("finalize", &finalize);
210 
211  declareMultiIndex(mod);
212  declareIndex(mod);
213  declareSolver(mod);
214 
215  return mod.ptr();
216 }
217 }
218 }
219 }
220 } // namespace lsst::meas::extensions::astrometryNet
lsst::afw::geom::Angle healpixDistance(int hp, int nside, lsst::afw::coord::Coord const &coord)
Calculate the distance from coordinates to a healpix.