SourceXtractorPlusPlus  0.14
Please provide a description of the project.
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
FitsImageSource.cpp
Go to the documentation of this file.
1 
17 /*
18  * FitsImageSource.cpp
19  *
20  * Created on: Mar 21, 2018
21  * Author: mschefer
22  */
23 
24 #include <iomanip>
25 #include <fstream>
26 #include <numeric>
27 #include <string>
28 
29 #include <boost/filesystem/path.hpp>
30 #include <boost/filesystem/operations.hpp>
31 #include <boost/regex.hpp>
32 #include <boost/algorithm/string/case_conv.hpp>
33 #include <boost/algorithm/string/trim.hpp>
34 
36 
38 #include "SEUtils/VariantCast.h"
39 
41 
42 namespace SourceXtractor {
43 
44 namespace {
45 }
46 
48  : m_filename(filename), m_manager(manager), m_hdu_number(hdu_number) {
49  int status = 0;
50  int bitpix, naxis;
51  long naxes[2] = {1,1};
52 
53  m_fits_file = m_manager->getFitsFile(filename);
54  auto fptr = m_fits_file->getFitsFilePtr();
55 
56  if (m_hdu_number <= 0) {
57  if (fits_get_hdu_num(fptr, &m_hdu_number) < 0) {
58  throw Elements::Exception() << "Can't get the active HDU from the FITS file: " << filename;
59  }
60  }
61  else {
62  switchHdu(fptr, m_hdu_number);
63  }
64 
65  fits_get_img_param(fptr, 2, &bitpix, &naxis, naxes, &status);
66  if (status != 0 || naxis != 2) {
67  throw Elements::Exception() << "Can't find 2D image in FITS file: " << filename << "[" << m_hdu_number << "]";
68  }
69 
70  m_width = naxes[0];
71  m_height = naxes[1];
72 
73  if (image_type < 0) {
74  switch (bitpix) {
75  case FLOAT_IMG:
77  break;
78  case DOUBLE_IMG:
80  break;
81  case LONG_IMG:
83  break;
84  case ULONG_IMG:
86  break;
87  case LONGLONG_IMG:
89  break;
90  default:
91  throw Elements::Exception() << "Unsupported FITS image type: " << bitpix;
92  }
93  } else {
94  m_image_type = image_type;
95  }
96 }
97 
98 
100  const std::shared_ptr<CoordinateSystem> coord_system, bool append, bool empty_primary,
102  : m_filename(filename), m_manager(manager), m_width(width), m_height(height), m_image_type(image_type) {
103  {
104 
105  m_manager->closeAndPurgeFile(filename);
106 
107  int status = 0;
108  fitsfile* fptr = nullptr;
109 
110  if (append) {
111  fits_open_image(&fptr, filename.c_str(), READWRITE, &status);
112  if (status != 0) {
113  status = 0;
114  fits_create_file(&fptr, filename.c_str(), &status);
115  if (empty_primary) {
116  fits_create_img(fptr, FLOAT_IMG, 0, nullptr, &status);
117  }
118  if (status != 0) {
119  throw Elements::Exception() << "Can't open or create FITS file to add HDU: " << filename;
120  }
121  }
122  } else {
123  // Overwrite file
124  fits_create_file(&fptr, ("!"+filename).c_str(), &status);
125  if (empty_primary) {
126  fits_create_img(fptr, FLOAT_IMG, 0, nullptr, &status);
127  }
128  if (status != 0) {
129  throw Elements::Exception() << "Can't create or overwrite FITS file: " << filename;
130  }
131  }
132  assert(fptr != nullptr);
133 
134  long naxes[2] = {width, height};
135  fits_create_img(fptr, getImageType(), 2, naxes, &status);
136 
137  if (fits_get_hdu_num(fptr, &m_hdu_number) < 0) {
138  throw Elements::Exception() << "Can't get the active HDU from the FITS file: " << filename;
139  }
140 
141  int hdutype = 0;
142  fits_movabs_hdu(fptr, m_hdu_number, &hdutype, &status);
143 
144  if (coord_system) {
145  auto headers = coord_system->getFitsHeaders();
146  for (auto& h : headers) {
147  std::ostringstream padded_key, serializer;
148  padded_key << std::setw(8) << std::left << h.first;
149 
150  serializer << padded_key.str() << "= " << std::left << std::setw(70) << h.second;
151  auto str = serializer.str();
152 
153  fits_update_card(fptr, padded_key.str().c_str(), str.c_str(), &status);
154  if (status != 0) {
155  char err_txt[31];
156  fits_get_errstatus(status, err_txt);
157  throw Elements::Exception() << "Couldn't write the WCS headers (" << err_txt << "): " << str;
158  }
159  }
160  }
161 
162  std::vector<char> buffer(width * ImageTile::getTypeSize(image_type));
163  for (int i = 0; i<height; i++) {
164  long first_pixel[2] = {1, i+1};
165  fits_write_pix(fptr, getDataType(), first_pixel, width, &buffer[0], &status);
166  }
167  fits_close_file(fptr, &status);
168 
169  if (status != 0) {
170  throw Elements::Exception() << "Couldn't allocate space for new FITS file: " << filename;
171  }
172  }
173 
174  // Reopens the newly created file through theFitsFileManager
175  m_fits_file = m_manager->getFitsFile(filename, true);
176  auto fptr = m_fits_file->getFitsFilePtr();
177  switchHdu(fptr, m_hdu_number);
178 }
179 
181  auto fptr = m_fits_file->getFitsFilePtr();
182  switchHdu(fptr, m_hdu_number);
183 
184  auto tile = ImageTile::create(m_image_type, x, y, width, height,
185  std::const_pointer_cast<ImageSource>(shared_from_this()));
186 
187  long first_pixel[2] = {x + 1, y + 1};
188  long last_pixel[2] = {x + width, y + height};
189  long increment[2] = {1, 1};
190  int status = 0;
191 
192  fits_read_subset(fptr, getDataType(), first_pixel, last_pixel, increment,
193  nullptr, tile->getDataPtr(), nullptr, &status);
194  if (status != 0) {
195  throw Elements::Exception() << "Error reading image tile from FITS file.";
196  }
197 
198  return tile;
199 }
200 
202  auto fptr = m_fits_file->getFitsFilePtr();
203  switchHdu(fptr, m_hdu_number);
204 
205  int x = tile.getPosX();
206  int y = tile.getPosY();
207  int width = tile.getWidth();
208  int height = tile.getHeight();
209 
210  long first_pixel[2] = {x + 1, y + 1};
211  long last_pixel[2] = {x + width, y + height};
212  int status = 0;
213 
214  fits_write_subset(fptr, getDataType(), first_pixel, last_pixel, tile.getDataPtr(), &status);
215  if (status != 0) {
216  throw Elements::Exception() << "Error saving image tile to FITS file.";
217  }
218 }
219 
220 void FitsImageSource::switchHdu(fitsfile *fptr, int hdu_number) const {
221  int status = 0;
222  int hdu_type = 0;
223 
224  fits_movabs_hdu(fptr, hdu_number, &hdu_type, &status);
225 
226  if (status != 0) {
227  throw Elements::Exception() << "Could not switch to HDU # " << hdu_number << " in file " << m_filename;
228  }
229  if (hdu_type != IMAGE_HDU) {
230  throw Elements::Exception() << "Trying to access non-image HDU in file " << m_filename;
231  }
232 }
233 
234 
235 
237  number_of_records = 0;
238  std::string records;
239 
240  auto& headers = getMetadata();
241  for (auto record : headers) {
242  auto key = record.first;
243 
244  std::string record_string(key);
245  if (record_string.size() > 8) {
246  throw Elements::Exception() << "FITS keyword longer than 8 characters";
247  } else if (record_string.size() < 8) {
248  record_string += std::string(8 - record_string.size(), ' ');
249  }
250 
251  if (headers.at(key).m_value.type() == typeid(std::string)) {
252  record_string += "= '" + VariantCast<std::string>(headers.at(key).m_value) + "'";
253  }
254  else {
255  record_string += "= " + VariantCast<std::string>(headers.at(key).m_value);
256  }
257 
258  if (record_string.size() > 80) {
259  throw Elements::Exception() << "FITS record longer than 80 characters";
260  }
261 
262 
263  if (record_string.size() < 80) {
264  record_string += std::string(80 - record_string.size(), ' ');
265  }
266 
267  records += record_string;
268  number_of_records++;
269  }
270 
271  std::string record_string("END");
272  record_string += std::string(80 - record_string.size(), ' ');
273  records += record_string;
274 
275  std::unique_ptr<std::vector<char>> buffer(new std::vector<char>(records.begin(), records.end()));
276  buffer->emplace_back(0);
277  return buffer;
278 }
279 
281  auto fptr = m_fits_file->getFitsFilePtr();
282  switchHdu(fptr, m_hdu_number);
283 
284  std::ostringstream padded_key, serializer;
285  padded_key << std::setw(8) << std::left << key;
286 
287  serializer << padded_key.str() << "= " << std::left << std::setw(70) << VariantCast<std::string>(value.m_value);
288  auto str = serializer.str();
289 
290  int status = 0;
291  fits_update_card(fptr, padded_key.str().c_str(), str.c_str(), &status);
292 
293  if (status != 0) {
294  char err_txt[31];
295  fits_get_errstatus(status, err_txt);
296  throw Elements::Exception() << "Couldn't write the metadata (" << err_txt << "): " << str;
297  }
298 
299  // update the metadata
300  m_fits_file->getHDUHeaders(m_hdu_number)[key] = value;
301 }
302 
304  switch (m_image_type) {
305  default:
307  return TFLOAT;
309  return TDOUBLE;
310  case ImageTile::IntImage:
311  return TINT;
313  return TUINT;
315  return TLONGLONG;
316  }
317 }
318 
320  switch (m_image_type) {
321  default:
323  return FLOAT_IMG;
325  return DOUBLE_IMG;
326  case ImageTile::IntImage:
327  return LONG_IMG;
329  return ULONG_IMG;
331  return LONGLONG_IMG;
332  }
333 }
334 
335 //FIXME add missing types
336 
337 }
std::unique_ptr< std::vector< char > > getFitsHeaders(int &number_of_records) const
const std::map< std::string, MetadataEntry > getMetadata() const override
std::shared_ptr< DependentParameter< std::shared_ptr< EngineParameter > > > x
T left(T...args)
ImageTile::ImageType m_image_type
T end(T...args)
std::shared_ptr< ImageTile > getImageTile(int x, int y, int width, int height) const override
std::shared_ptr< DependentParameter< std::shared_ptr< EngineParameter > > > y
std::shared_ptr< FitsFile > m_fits_file
T setw(T...args)
int getPosX() const
Definition: ImageTile.h:56
STL class.
std::shared_ptr< FitsFileManager > m_manager
string filename
Definition: conf.py:63
static size_t getTypeSize(ImageType image_type)
Definition: ImageTile.h:119
void setMetadata(std::string key, MetadataEntry value) override
int getWidth() const
Definition: ImageTile.h:66
void switchHdu(fitsfile *fptr, int hdu_number) const
T size(T...args)
STL class.
STL class.
static std::shared_ptr< ImageTile > create(ImageType image_type, int x, int y, int width, int height, std::shared_ptr< ImageSource > source=nullptr)
Definition: ImageTile.cpp:24
T begin(T...args)
int getHeight() const
Definition: ImageTile.h:70
T c_str(T...args)
void saveTile(ImageTile &tile) override
int getPosY() const
Definition: ImageTile.h:60
virtual void * getDataPtr()=0
FitsImageSource(const std::string &filename, int hdu_number=0, ImageTile::ImageType image_type=ImageTile::AutoType, std::shared_ptr< FitsFileManager > manager=FitsFileManager::getInstance())