Alexandria 2.31.0
SDC-CH common library for the Euclid project
Loading...
Searching...
No Matches
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
35using Euclid::Table::operator<<;
36
37namespace Euclid {
38namespace Table {
39
40using NdArray::NdArray;
41
42template <typename T>
44 std::ostringstream stream;
45 stream << std::scientific << cell_stream_adaptor(value);
46 return stream.str();
47}
48
49size_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
57size_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
88template <typename T>
89size_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
102template <typename T>
103size_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
121template <typename T>
122static 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
131template <typename T>
132static 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
137template <typename T>
138static 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
188template <typename T>
189std::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
196template <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
206template <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
216template <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
228template <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
241template <typename T>
242void 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
251template <typename T>
252void 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
261std::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;
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
300void 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
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< T > createSingleNdArrayVectorColumnData(const Euclid::Table::Table &table, size_t column_index)
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)
std::vector< T > createSingleValueVectorColumnData(const Euclid::Table::Table &table, size_t column_index)
void populateColumn(const Table &table, int column_index, const CCfits::ExtHDU &table_hdu, long first_row)
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< std::valarray< T > > createNdArrayColumnData(const Euclid::Table::Table &table, size_t column_index)
std::string scientificFormat(const T &value)
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::vector< T > createColumnData(const Table &table, size_t column_index)
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)