SourceXtractorPlusPlus  0.11
Please provide a description of the project.
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
WCS.cpp
Go to the documentation of this file.
1 
17 /*
18  * WCS.cpp
19  *
20  * Created on: Nov 17, 2016
21  * Author: mschefer
22  */
23 
24 #include <mutex>
25 
26 #include <boost/algorithm/string/trim.hpp>
27 
28 #include <fitsio.h>
29 
30 #include <wcslib/wcs.h>
31 #include <wcslib/wcshdr.h>
32 #include <wcslib/wcsfix.h>
33 #include <wcslib/wcsprintf.h>
34 #include <wcslib/getwcstab.h>
35 
37 #include "ElementsKernel/Logging.h"
38 
40 
41 namespace SourceXtractor {
42 
43 static auto logger = Elements::Logging::getLogger("WCS");
44 
45 decltype(&lincpy) safe_lincpy = &lincpy;
46 
52 static int wrapped_lincpy(int alloc, const struct linprm *linsrc, struct linprm *lindst) {
53  static std::mutex lincpy_mutex;
54  std::lock_guard<std::mutex> lock(lincpy_mutex);
55 
56  return lincpy(alloc, linsrc, lindst);
57 }
58 
59 
60 WCS::WCS(const FitsImageSource<SeFloat>& fits_image_source) : m_wcs(nullptr, nullptr) {
61  int number_of_records = 0;
62  auto fits_headers = fits_image_source.getFitsHeaders(number_of_records);
63 
64  int nreject = 0, nwcs = 0;
65  wcsprm* wcs;
66  wcspih(&(*fits_headers)[0], number_of_records, WCSHDR_all, 0, &nreject, &nwcs, &wcs);
67  wcsset(wcs);
68 
69  m_wcs = decltype(m_wcs)(wcs, [nwcs](wcsprm* wcs) {
70  int nwcs_copy = nwcs;
71  wcsfree(wcs);
72  wcsvfree(&nwcs_copy, &wcs);
73  });
74 
75 
77 }
78 
80 }
81 
83  // wcsprm is in/out, since its member lin is modified by wcsp2s
84  wcsprm wcs_copy = *m_wcs;
85  wcs_copy.lin.flag = -1;
86  safe_lincpy(true, &m_wcs->lin, &wcs_copy.lin);
87  linset(&wcs_copy.lin);
88 
89  // +1 as fits standard coordinates start at 1
90  double pc_array[2] {image_coordinate.m_x + 1, image_coordinate.m_y + 1};
91 
92  double ic_array[2] {0, 0};
93  double wc_array[2] {0, 0};
94  double phi, theta;
95 
96  int status = 0;
97  wcsp2s(&wcs_copy, 1, 1, pc_array, ic_array, &phi, &theta, wc_array, &status);
98  int ret_val = wcsp2s(&wcs_copy, 1, 1, pc_array, ic_array, &phi, &theta, wc_array, &status);
99  linfree(&wcs_copy.lin);
100  if (ret_val != 0) {
101  logger.error() << "wcslib's wcsp2s returned with error code: " << ret_val;
102  throw Elements::Exception() << "WCS exception";
103  }
104 
105  return WorldCoordinate(wc_array[0], wc_array[1]);
106 }
107 
109  // wcsprm is in/out, since its member lin is modified by wcss2p
110  wcsprm wcs_copy = *m_wcs;
111  wcs_copy.lin.flag = -1;
112  safe_lincpy(true, &m_wcs->lin, &wcs_copy.lin);
113  linset(&wcs_copy.lin);
114 
115  double pc_array[2] {0, 0};
116  double ic_array[2] {0, 0};
117  double wc_array[2] {world_coordinate.m_alpha, world_coordinate.m_delta};
118  double phi, theta;
119 
120  int status = 0;
121  int ret_val = wcss2p(&wcs_copy, 1, 1, wc_array, &phi, &theta, ic_array, pc_array, &status);
122  linfree(&wcs_copy.lin);
123  if (ret_val != 0) {
124  logger.error() << "wcslib's wcss2p returned with error code: " << ret_val;
125  throw Elements::Exception() << "WCS exception";
126  }
127 
128  return ImageCoordinate(pc_array[0] - 1, pc_array[1] - 1); // -1 as fits standard coordinates start at 1
129 }
130 
132  int nkeyrec;
133  char *raw_header;
134 
135  if (wcshdo(WCSHDO_none, m_wcs.get(), &nkeyrec, &raw_header) != 0) {
136  throw Elements::Exception() << "Failed to get the FITS headers for the WCS coordinate system";
137  }
138 
140  for (int i = 0; i < nkeyrec; ++i) {
141  char *hptr = &raw_header[80 * i];
142  std::string key(hptr, hptr + 8);
143  boost::trim(key);
144  std::string value(hptr + 9, hptr + 72);
145  boost::trim(value);
146  if (!key.empty()) {
147  headers.emplace(std::make_pair(key, value));
148  }
149  }
150 
151  free(raw_header);
152  return headers;
153 }
154 
155 }
std::unique_ptr< std::vector< char > > getFitsHeaders(int &number_of_records) const
T empty(T...args)
static Elements::Logging logger
WCS(const FitsImageSource< SeFloat > &fits_image_source)
Definition: WCS.cpp:60
decltype(&lincpy) safe_lincpy
Definition: WCS.cpp:45
T free(T...args)
STL class.
WorldCoordinate imageToWorld(ImageCoordinate image_coordinate) const override
Definition: WCS.cpp:82
T lock(T...args)
T make_pair(T...args)
std::map< std::string, std::string > getFitsHeaders() const override
Definition: WCS.cpp:131
std::unique_ptr< wcsprm, std::function< void(wcsprm *)> > m_wcs
Definition: WCS.h:48
virtual ~WCS()
Definition: WCS.cpp:79
static int wrapped_lincpy(int alloc, const struct linprm *linsrc, struct linprm *lindst)
Definition: WCS.cpp:52
void error(const std::string &logMessage)
T emplace(T...args)
ImageCoordinate worldToImage(WorldCoordinate world_coordinate) const override
Definition: WCS.cpp:108
static Logging getLogger(const std::string &name="")