Alexandria  2.27.0
SDC-CH common library for the Euclid project
FitsWriterHelper.cpp
Go to the documentation of this file.
1 
25 #include "FitsWriterHelper.h"
27 #include "ReaderHelper.h"
28 #include "Table/Table.h"
29 #include <CCfits/CCfits>
30 #include <algorithm>
31 #include <boost/lexical_cast.hpp>
32 #include <sstream>
33 #include <valarray>
34 
35 using Euclid::Table::operator<<;
36 
37 namespace Euclid {
38 namespace Table {
39 
40 using NdArray::NdArray;
41 
42 template <typename T>
43 std::string scientificFormat(const T& value) {
44  std::ostringstream stream;
45  stream << std::scientific << cell_stream_adaptor(value);
46  return stream.str();
47 }
48 
49 size_t maxWidth(const Table& table, size_t column_index) {
50  size_t width = table.getColumnInfo()->getDescription(column_index).size;
51  for (const auto& row : table) {
52  width = std::max(width, boost::lexical_cast<std::string>(cell_stream_adaptor(row[column_index])).size());
53  }
54  return width;
55 }
56 
57 size_t maxWidthScientific(const Table& table, size_t column_index) {
58  size_t width = 0;
59  for (const auto& row : table) {
60  width = std::max(width, scientificFormat(row[column_index]).size());
61  }
62  return width;
63 }
64 
66  auto column_info = table.getColumnInfo();
67  std::vector<std::string> format_list{};
68  for (size_t column_index = 0; column_index < column_info->size(); ++column_index) {
69  auto type = column_info->getDescription(column_index).type;
70  if (type == typeid(bool)) {
71  format_list.emplace_back("I1");
72  } else if (type == typeid(int32_t) || type == typeid(int64_t)) {
73  size_t width = maxWidth(table, column_index);
74  format_list.emplace_back("I" + boost::lexical_cast<std::string>(std::max(static_cast<size_t>(1), width)));
75  } else if (type == typeid(float) || type == typeid(double)) {
76  size_t width = maxWidthScientific(table, column_index);
77  format_list.emplace_back("E" + boost::lexical_cast<std::string>(std::max(static_cast<size_t>(1), width)));
78  } else if (type == typeid(std::string)) {
79  size_t width = maxWidth(table, column_index);
80  format_list.emplace_back("A" + boost::lexical_cast<std::string>(std::max(static_cast<size_t>(1), width)));
81  } else {
82  throw Elements::Exception() << "Unsupported column format for FITS ASCII table export: " << type.name();
83  }
84  }
85  return format_list;
86 }
87 
88 template <typename T>
89 size_t vectorSize(const Table& table, size_t column_index) {
90  if (table.size() == 0) {
91  return 0;
92  }
93  size_t size = boost::get<std::vector<T>>(table[0][column_index]).size();
94  for (const auto& row : table) {
95  if (boost::get<std::vector<T>>(row[column_index]).size() != size) {
96  throw Elements::Exception() << "Binary FITS table variable length vector columns are not supported";
97  }
98  }
99  return size;
100 }
101 
102 template <typename T>
103 size_t ndArraySize(const Table& table, size_t column_index) {
104  if (table.size() == 0) {
105  return 0;
106  }
107  const auto& ndarray = boost::get<NdArray<T>>(table[0][column_index]);
108  size_t size = ndarray.size();
109  auto shape = ndarray.shape();
110  for (const auto& row : table) {
111  if (boost::get<NdArray<T>>(row[column_index]).shape() != shape) {
112  throw Elements::Exception() << "Binary FITS table variable shape array columns are not supported";
113  }
114  }
115  return size;
116 }
117 
118 // Defined in FitsReaderHelper.cpp
120 
121 template <typename T>
122 static std::string GenericScalarFormat(const Table&, size_t) {
123  auto i = std::find_if(ScalarTypeMap.begin(), ScalarTypeMap.end(),
124  [](const std::pair<char, std::type_index>& p) { return p.second == typeid(T); });
125  if (i == ScalarTypeMap.end()) {
126  throw Elements::Exception() << "Unsupported column format for FITS binary table export: " << typeid(T).name();
127  }
128  return std::string(1, i->first);
129 }
130 
131 template <typename T>
132 static std::string GenericVectorFormat(const Table& table, size_t column_index) {
133  size_t size = vectorSize<T>(table, column_index);
134  return boost::lexical_cast<std::string>(size) + GenericScalarFormat<T>(table, column_index);
135 }
136 
137 template <typename T>
138 static std::string GenericNdFormat(const Table& table, size_t column_index) {
139  size_t size = ndArraySize<T>(table, column_index);
140  return boost::lexical_cast<std::string>(size) + GenericScalarFormat<T>(table, column_index);
141 }
142 
144  // Scalars
145  {typeid(bool), GenericScalarFormat<bool>},
146  {typeid(int32_t), GenericScalarFormat<int32_t>},
147  {typeid(int64_t), GenericScalarFormat<int64_t>},
148  {typeid(float), GenericScalarFormat<float>},
149  {typeid(double), GenericScalarFormat<double>},
150  // Vectors
151  {typeid(std::vector<bool>), GenericVectorFormat<bool>},
152  {typeid(std::vector<int32_t>), GenericVectorFormat<int32_t>},
153  {typeid(std::vector<int64_t>), GenericVectorFormat<int64_t>},
154  {typeid(std::vector<float>), GenericVectorFormat<float>},
155  {typeid(std::vector<double>), GenericVectorFormat<double>},
156  // String
157  {typeid(std::string),
158  [](const Table& table, size_t column) {
159  size_t width = maxWidth(table, column);
160  // CCfits uses the width as denominator, so it can't be 0!
161  return boost::lexical_cast<std::string>(std::max(static_cast<size_t>(1), width)) + "A";
162  }},
163  // NdArray
164  {typeid(NdArray<int32_t>), GenericNdFormat<int32_t>},
165  {typeid(NdArray<int64_t>), GenericNdFormat<int64_t>},
166  {typeid(NdArray<float>), GenericNdFormat<float>},
167  {typeid(NdArray<double>), GenericNdFormat<double>},
168 };
169 
171  auto column_info = table.getColumnInfo();
172  std::vector<std::string> format_list;
173  format_list.reserve(column_info->size());
174 
175  for (size_t column_index = 0; column_index < column_info->size(); ++column_index) {
176  auto type = column_info->getDescription(column_index).type;
177 
178  auto i = BinaryFormatter.find(type);
179  if (i == BinaryFormatter.end()) {
180  throw Elements::Exception() << "Unsupported column format for FITS binary table export: " << type.name();
181  }
182 
183  format_list.emplace_back(i->second(table, column_index));
184  }
185  return format_list;
186 }
187 
188 template <typename T>
189 std::vector<T> createColumnData(const Table& table, size_t column_index) {
190  std::vector<T> data(table.size());
191  std::transform(table.begin(), table.end(), data.begin(),
192  [column_index](const Row& row) { return boost::get<T>(row[column_index]); });
193  return data;
194 }
195 
196 template <typename T>
199  for (auto& row : table) {
200  const auto& vec = boost::get<std::vector<T>>(row[column_index]);
201  result.emplace_back(vec.data(), vec.size());
202  }
203  return result;
204 }
205 
206 template <typename T>
208  std::vector<T> result{};
209  for (auto& row : table) {
210  const auto& vec = boost::get<std::vector<T>>(row[column_index]);
211  result.push_back(vec.front());
212  }
213  return result;
214 }
215 
216 template <typename T>
219  for (auto& row : table) {
220  const auto& ndarray = boost::get<NdArray<T>>(row[column_index]);
221  std::valarray<T> data(ndarray.size());
222  std::copy(std::begin(ndarray), std::end(ndarray), std::begin(data));
223  result.emplace_back(std::move(data));
224  }
225  return result;
226 }
227 
228 template <typename T>
230  std::vector<T> result{};
231  for (auto& row : table) {
232  const auto& nd = boost::get<NdArray<T>>(row[column_index]);
233  if (nd.size() > 0)
234  result.push_back(*nd.begin());
235  else
236  result.push_back(0);
237  }
238  return result;
239 }
240 
241 template <typename T>
242 void populateVectorColumn(const Table& table, int column_index, const CCfits::ExtHDU& table_hdu, long first_row) {
243  const auto& vec = boost::get<std::vector<T>>(table[0][column_index]);
244  if (vec.size() > 1) {
245  table_hdu.column(column_index + 1).writeArrays(createVectorColumnData<T>(table, column_index), first_row);
246  } else {
247  table_hdu.column(column_index + 1).write(createSingleValueVectorColumnData<T>(table, column_index), first_row);
248  }
249 }
250 
251 template <typename T>
252 void populateNdArrayColumn(const Table& table, int column_index, const CCfits::ExtHDU& table_hdu, long first_row) {
253  const auto& nd = boost::get<NdArray<T>>(table[0][column_index]);
254  if (nd.size() > 1) {
255  table_hdu.column(column_index + 1).writeArrays(createNdArrayColumnData<T>(table, column_index), first_row);
256  } else {
257  table_hdu.column(column_index + 1).write(createSingleNdArrayVectorColumnData<T>(table, column_index), first_row);
258  }
259 }
260 
261 std::string getTDIM(const Table& table, int column_index) {
262  if (table.size() == 0) {
263  return "";
264  }
265 
266  auto& first_row = table[0];
267  auto& cell = first_row[column_index];
268  auto type = table.getColumnInfo()->getDescription(column_index).type;
269  std::vector<size_t> shape;
270 
271  if (type == typeid(NdArray<int32_t>)) {
272  shape = boost::get<NdArray<int32_t>>(cell).shape();
273  } else if (type == typeid(NdArray<int64_t>)) {
274  shape = boost::get<NdArray<int64_t>>(cell).shape();
275  } else if (type == typeid(NdArray<float>)) {
276  shape = boost::get<NdArray<float>>(cell).shape();
277  } else if (type == typeid(NdArray<double>)) {
278  shape = boost::get<NdArray<double>>(cell).shape();
279  } else {
280  return "";
281  }
282 
283  int64_t ncells = std::accumulate(shape.begin(), shape.end(), 1, std::multiplies<size_t>());
284  if (ncells == 1) {
285  return "";
286  }
287 
288  std::stringstream stream;
289  stream << '(';
290 
291  size_t j;
292  for (j = shape.size() - 1; j > 0; --j) {
293  stream << shape[j] << ",";
294  }
295 
296  stream << shape[j] << ')';
297  return stream.str();
298 }
299 
300 void populateColumn(const Table& table, int column_index, const CCfits::ExtHDU& table_hdu, long first_row) {
301  if (table.size() == 0) {
302  return;
303  }
304  auto type = table.getColumnInfo()->getDescription(column_index).type;
305  // CCfits indices start from 1
306  if (type == typeid(bool)) {
307  table_hdu.column(column_index + 1).write(createColumnData<bool>(table, column_index), first_row);
308  } else if (type == typeid(int32_t)) {
309  table_hdu.column(column_index + 1).write(createColumnData<int32_t>(table, column_index), first_row);
310  } else if (type == typeid(int64_t)) {
311  table_hdu.column(column_index + 1).write(createColumnData<int64_t>(table, column_index), first_row);
312  } else if (type == typeid(float)) {
313  table_hdu.column(column_index + 1).write(createColumnData<float>(table, column_index), first_row);
314  } else if (type == typeid(double)) {
315  table_hdu.column(column_index + 1).write(createColumnData<double>(table, column_index), first_row);
316  } else if (type == typeid(std::string)) {
317  table_hdu.column(column_index + 1).write(createColumnData<std::string>(table, column_index), first_row);
318  } else if (type == typeid(std::vector<int32_t>)) {
319  populateVectorColumn<int32_t>(table, column_index, table_hdu, first_row);
320  } else if (type == typeid(std::vector<int64_t>)) {
321  populateVectorColumn<int64_t>(table, column_index, table_hdu, first_row);
322  } else if (type == typeid(std::vector<float>)) {
323  populateVectorColumn<float>(table, column_index, table_hdu, first_row);
324  } else if (type == typeid(std::vector<double>)) {
325  populateVectorColumn<double>(table, column_index, table_hdu, first_row);
326  } else if (type == typeid(NdArray<int32_t>)) {
327  populateNdArrayColumn<int32_t>(table, column_index, table_hdu, first_row);
328  } else if (type == typeid(NdArray<int64_t>)) {
329  populateNdArrayColumn<int64_t>(table, column_index, table_hdu, first_row);
330  } else if (type == typeid(NdArray<float>)) {
331  populateNdArrayColumn<float>(table, column_index, table_hdu, first_row);
332  } else if (type == typeid(NdArray<double>)) {
333  populateNdArrayColumn<double>(table, column_index, table_hdu, first_row);
334  } else {
335  throw Elements::Exception() << "Cannot populate FITS column with data of type " << type.name();
336  }
337 }
338 
339 } // namespace Table
340 } // end of namespace Euclid
T accumulate(T... args)
T begin(T... args)
const std::vector< size_t > & shape() const
Definition: NdArray.h:319
NdArray(std::vector< size_t > shape_)
Represents one row of a Table.
Definition: Row.h:57
Represents a table.
Definition: Table.h:49
T copy(T... args)
T emplace_back(T... args)
T end(T... args)
T find_if(T... args)
T scientific(T... args)
T max(T... args)
T move(T... args)
const std::map< std::type_index, std::function< std::string(const Table &, size_t)> > BinaryFormatter
std::vector< std::string > getAsciiFormatList(const Table &table)
Returns a vector with strings representing the FITS ASCII table formats for the given table.
std::vector< std::valarray< T > > createVectorColumnData(const Euclid::Table::Table &table, size_t column_index)
static std::string GenericVectorFormat(const Table &table, size_t column_index)
void populateColumn(const Table &table, int column_index, const CCfits::ExtHDU &table_hdu, long first_row)
std::vector< T > createSingleValueVectorColumnData(const Euclid::Table::Table &table, size_t column_index)
static std::string GenericNdFormat(const Table &table, size_t column_index)
size_t maxWidthScientific(const Table &table, size_t column_index)
size_t maxWidth(const Table &table, size_t column_index)
const std::vector< std::pair< char, std::type_index > > ScalarTypeMap
size_t ndArraySize(const Table &table, size_t column_index)
std::vector< std::string > getBinaryFormatList(const Table &table)
Returns a vector with strings representing the FITS binary table formats for the given table.
std::vector< T > createColumnData(const Table &table, size_t column_index)
std::vector< T > createSingleNdArrayVectorColumnData(const Euclid::Table::Table &table, size_t column_index)
std::string scientificFormat(const T &value)
std::vector< std::valarray< T > > createNdArrayColumnData(const Euclid::Table::Table &table, size_t column_index)
size_t vectorSize(const Table &table, size_t column_index)
void populateVectorColumn(const Table &table, int column_index, const CCfits::ExtHDU &table_hdu, long first_row)
static std::string GenericScalarFormat(const Table &, size_t)
std::string getTDIM(const Table &table, int column_index)
void populateNdArrayColumn(const Table &table, int column_index, const CCfits::ExtHDU &table_hdu, long first_row)
T push_back(T... args)
T reserve(T... args)
T size(T... args)
T str(T... args)
T transform(T... args)