Alexandria 2.25.0
SDC-CH common library for the Euclid project
FitsWriterHelper.cpp
Go to the documentation of this file.
1/*
2 * Copyright (C) 2012-2022 Euclid Science Ground Segment
3 *
4 * This library is free software; you can redistribute it and/or modify it under
5 * the terms of the GNU Lesser General Public License as published by the Free
6 * Software Foundation; either version 3.0 of the License, or (at your option)
7 * any later version.
8 *
9 * This library is distributed in the hope that it will be useful, but WITHOUT
10 * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
11 * FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License for more
12 * details.
13 *
14 * You should have received a copy of the GNU Lesser General Public License
15 * along with this library; if not, write to the Free Software Foundation, Inc.,
16 * 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA
17 */
18
25#include "FitsWriterHelper.h"
27#include "Table/Table.h"
28#include <CCfits/CCfits>
29#include <algorithm>
30#include <boost/lexical_cast.hpp>
31#include <iomanip>
32#include <sstream>
33#include <valarray>
34
35namespace Euclid {
36namespace Table {
37
39
40template <typename T>
42 std::ostringstream stream;
43 stream << std::scientific << value;
44 return stream.str();
45}
46
47size_t maxWidth(const Table& table, size_t column_index) {
48 size_t width = table.getColumnInfo()->getDescription(column_index).size;
49 for (const auto& row : table) {
50 width = std::max(width, boost::lexical_cast<std::string>(row[column_index]).size());
51 }
52 return width;
53}
54
55size_t maxWidthScientific(const Table& table, size_t column_index) {
56 size_t width = 0;
57 for (const auto& row : table) {
58 width = std::max(width, scientificFormat(row[column_index]).size());
59 }
60 return width;
61}
62
64 auto column_info = table.getColumnInfo();
65 std::vector<std::string> format_list{};
66 for (size_t column_index = 0; column_index < column_info->size(); ++column_index) {
67 auto type = column_info->getDescription(column_index).type;
68 if (type == typeid(bool)) {
69 format_list.push_back("I1");
70 } else if (type == typeid(int32_t) || type == typeid(int64_t)) {
71 size_t width = maxWidth(table, column_index);
72 format_list.push_back("I" + boost::lexical_cast<std::string>(std::max(static_cast<size_t>(1), width)));
73 } else if (type == typeid(float) || type == typeid(double)) {
74 size_t width = maxWidthScientific(table, column_index);
75 format_list.push_back("E" + boost::lexical_cast<std::string>(std::max(static_cast<size_t>(1), width)));
76 } else if (type == typeid(std::string)) {
77 size_t width = maxWidth(table, column_index);
78 format_list.push_back("A" + boost::lexical_cast<std::string>(std::max(static_cast<size_t>(1), width)));
79 } else {
80 throw Elements::Exception() << "Unsupported column format for FITS ASCII table export: " << type.name();
81 }
82 }
83 return format_list;
84}
85
86template <typename T>
87size_t vectorSize(const Table& table, size_t column_index) {
88 if (table.size() == 0) {
89 return 0;
90 }
91 size_t size = boost::get<std::vector<T>>(table[0][column_index]).size();
92 for (const auto& row : table) {
93 if (boost::get<std::vector<T>>(row[column_index]).size() != size) {
94 throw Elements::Exception() << "Binary FITS table variable length vector columns are not supported";
95 }
96 }
97 return size;
98}
99
100template <typename T>
101size_t ndArraySize(const Table& table, size_t column_index) {
102 if (table.size() == 0) {
103 return 0;
104 }
105 const auto& ndarray = boost::get<NdArray<T>>(table[0][column_index]);
106 size_t size = ndarray.size();
107 auto shape = ndarray.shape();
108 for (const auto& row : table) {
109 if (boost::get<NdArray<T>>(row[column_index]).shape() != shape) {
110 throw Elements::Exception() << "Binary FITS table variable shape array columns are not supported";
111 }
112 }
113 return size;
114}
115
116// Defined in FitsReaderHelper.cpp
118
119template <typename T>
120static std::string GenericScalarFormat(const Table&, size_t) {
121 char type_format = '\0';
122 for (auto p = ScalarTypeMap.begin(); p != ScalarTypeMap.end(); ++p) {
123 if (p->second == typeid(T)) {
124 type_format = p->first;
125 break;
126 }
127 }
128 if (!type_format) {
129 throw Elements::Exception() << "Unsupported column format for FITS binary table export: " << typeid(T).name();
130 }
131 return std::string(1, type_format);
132}
133
134template <typename T>
135static std::string GenericVectorFormat(const Table& table, size_t column_index) {
136 size_t size = vectorSize<T>(table, column_index);
137 return boost::lexical_cast<std::string>(size) + GenericScalarFormat<T>(table, column_index);
138}
139
140template <typename T>
141static std::string GenericNdFormat(const Table& table, size_t column_index) {
142 size_t size = ndArraySize<T>(table, column_index);
143 return boost::lexical_cast<std::string>(size) + GenericScalarFormat<T>(table, column_index);
144}
145
147 // Scalars
148 {typeid(bool), GenericScalarFormat<bool>},
149 {typeid(int32_t), GenericScalarFormat<int32_t>},
150 {typeid(int64_t), GenericScalarFormat<int64_t>},
151 {typeid(float), GenericScalarFormat<float>},
152 {typeid(double), GenericScalarFormat<double>},
153 // Vectors
154 {typeid(std::vector<bool>), GenericVectorFormat<bool>},
155 {typeid(std::vector<int32_t>), GenericVectorFormat<int32_t>},
156 {typeid(std::vector<int64_t>), GenericVectorFormat<int64_t>},
157 {typeid(std::vector<float>), GenericVectorFormat<float>},
158 {typeid(std::vector<double>), GenericVectorFormat<double>},
159 // String
160 {typeid(std::string),
161 [](const Table& table, size_t column) {
162 size_t width = maxWidth(table, column);
163 // CCfits uses the width as denominator, so it can't be 0!
164 return boost::lexical_cast<std::string>(std::max(static_cast<size_t>(1), width)) + "A";
165 }},
166 // NdArray
167 {typeid(NdArray<int32_t>), GenericNdFormat<int32_t>},
168 {typeid(NdArray<int64_t>), GenericNdFormat<int64_t>},
169 {typeid(NdArray<float>), GenericNdFormat<float>},
170 {typeid(NdArray<double>), GenericNdFormat<double>},
171};
172
174 auto column_info = table.getColumnInfo();
175 std::vector<std::string> format_list;
176 format_list.reserve(column_info->size());
177
178 for (size_t column_index = 0; column_index < column_info->size(); ++column_index) {
179 auto type = column_info->getDescription(column_index).type;
180
181 auto i = BinaryFormatter.find(type);
182 if (i == BinaryFormatter.end()) {
183 throw Elements::Exception() << "Unsupported column format for FITS binary table export: " << type.name();
184 }
185
186 format_list.emplace_back(i->second(table, column_index));
187 }
188 return format_list;
189}
190
191template <typename T>
192std::vector<T> createColumnData(const Euclid::Table::Table& table, size_t column_index) {
193 std::vector<T> data{};
194 for (const auto& row : table) {
195 data.push_back(boost::get<T>(row[column_index]));
196 }
197 return data;
198}
199
200template <typename T>
203 for (auto& row : table) {
204 const auto& vec = boost::get<std::vector<T>>(row[column_index]);
205 result.emplace_back(vec.data(), vec.size());
206 }
207 return result;
208}
209
210template <typename T>
212 std::vector<T> result{};
213 for (auto& row : table) {
214 const auto& vec = boost::get<std::vector<T>>(row[column_index]);
215 result.push_back(vec.front());
216 }
217 return result;
218}
219
220template <typename T>
223 for (auto& row : table) {
224 const auto& ndarray = boost::get<NdArray<T>>(row[column_index]);
225 std::valarray<T> data(ndarray.size());
226 std::copy(std::begin(ndarray), std::end(ndarray), std::begin(data));
227 result.emplace_back(std::move(data));
228 }
229 return result;
230}
231
232template <typename T>
234 std::vector<T> result{};
235 for (auto& row : table) {
236 const auto& nd = boost::get<NdArray<T>>(row[column_index]);
237 if (nd.size() > 0)
238 result.push_back(*nd.begin());
239 else
240 result.push_back(0);
241 }
242 return result;
243}
244
245template <typename T>
246void populateVectorColumn(const Table& table, size_t column_index, CCfits::ExtHDU& table_hdu, long first_row) {
247 const auto& vec = boost::get<std::vector<T>>(table[0][column_index]);
248 if (vec.size() > 1) {
249 table_hdu.column(column_index + 1).writeArrays(createVectorColumnData<T>(table, column_index), first_row);
250 } else {
251 table_hdu.column(column_index + 1).write(createSingleValueVectorColumnData<T>(table, column_index), first_row);
252 }
253}
254
255template <typename T>
256void populateNdArrayColumn(const Table& table, size_t column_index, CCfits::ExtHDU& table_hdu, long first_row) {
257 const auto& nd = boost::get<NdArray<T>>(table[0][column_index]);
258 if (nd.size() > 1) {
259 table_hdu.column(column_index + 1).writeArrays(createNdArrayColumnData<T>(table, column_index), first_row);
260 } else {
261 table_hdu.column(column_index + 1).write(createSingleNdArrayVectorColumnData<T>(table, column_index), first_row);
262 }
263}
264
265std::string getTDIM(const Table& table, size_t column_index) {
266 if (table.size() == 0) {
267 return "";
268 }
269
270 auto& first_row = table[0];
271 auto& cell = first_row[column_index];
272 auto type = table.getColumnInfo()->getDescription(column_index).type;
274
275 if (type == typeid(NdArray<int32_t>)) {
276 shape = boost::get<NdArray<int32_t>>(cell).shape();
277 } else if (type == typeid(NdArray<int64_t>)) {
278 shape = boost::get<NdArray<int64_t>>(cell).shape();
279 } else if (type == typeid(NdArray<float>)) {
280 shape = boost::get<NdArray<float>>(cell).shape();
281 } else if (type == typeid(NdArray<double>)) {
282 shape = boost::get<NdArray<double>>(cell).shape();
283 } else {
284 return "";
285 }
286
287 int64_t ncells = 1;
288 for (auto& axis : shape) {
289 ncells *= axis;
290 }
291 if (ncells == 1) {
292 return "";
293 }
294
295 std::stringstream stream;
296 stream << '(';
297
298 int j;
299 for (j = shape.size() - 1; j > 0; --j) {
300 stream << shape[j] << ",";
301 }
302
303 stream << shape[j] << ')';
304 return stream.str();
305}
306
307void populateColumn(const Table& table, size_t column_index, CCfits::ExtHDU& table_hdu, long first_row) {
308 if (table.size() == 0) {
309 return;
310 }
311 auto type = table.getColumnInfo()->getDescription(column_index).type;
312 // CCfits indices start from 1
313 if (type == typeid(bool)) {
314 table_hdu.column(column_index + 1).write(createColumnData<bool>(table, column_index), first_row);
315 } else if (type == typeid(int32_t)) {
316 table_hdu.column(column_index + 1).write(createColumnData<int32_t>(table, column_index), first_row);
317 } else if (type == typeid(int64_t)) {
318 table_hdu.column(column_index + 1).write(createColumnData<int64_t>(table, column_index), first_row);
319 } else if (type == typeid(float)) {
320 table_hdu.column(column_index + 1).write(createColumnData<float>(table, column_index), first_row);
321 } else if (type == typeid(double)) {
322 table_hdu.column(column_index + 1).write(createColumnData<double>(table, column_index), first_row);
323 } else if (type == typeid(std::string)) {
324 table_hdu.column(column_index + 1).write(createColumnData<std::string>(table, column_index), first_row);
325 } else if (type == typeid(std::vector<int32_t>)) {
326 populateVectorColumn<int32_t>(table, column_index, table_hdu, first_row);
327 } else if (type == typeid(std::vector<int64_t>)) {
328 populateVectorColumn<int64_t>(table, column_index, table_hdu, first_row);
329 } else if (type == typeid(std::vector<float>)) {
330 populateVectorColumn<float>(table, column_index, table_hdu, first_row);
331 } else if (type == typeid(std::vector<double>)) {
332 populateVectorColumn<double>(table, column_index, table_hdu, first_row);
333 } else if (type == typeid(NdArray<int32_t>)) {
334 populateNdArrayColumn<int32_t>(table, column_index, table_hdu, first_row);
335 } else if (type == typeid(NdArray<int64_t>)) {
336 populateNdArrayColumn<int64_t>(table, column_index, table_hdu, first_row);
337 } else if (type == typeid(NdArray<float>)) {
338 populateNdArrayColumn<float>(table, column_index, table_hdu, first_row);
339 } else if (type == typeid(NdArray<double>)) {
340 populateNdArrayColumn<double>(table, column_index, table_hdu, first_row);
341 } else {
342 throw Elements::Exception() << "Cannot populate FITS column with data of type " << type.name();
343 }
344}
345
346} // namespace Table
347} // end of namespace Euclid
T begin(T... args)
const std::vector< size_t > & shape() const
Definition: NdArray.h:317
NdArray(std::vector< size_t > shape_)
Represents a table.
Definition: Table.h:49
T copy(T... args)
T emplace_back(T... args)
T end(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, size_t column_index, CCfits::ExtHDU &table_hdu, long first_row)
static std::string GenericNdFormat(const Table &table, size_t column_index)
std::string getTDIM(const Table &table, size_t column_index)
std::string scientificFormat(T value)
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.
void populateVectorColumn(const Table &table, size_t column_index, CCfits::ExtHDU &table_hdu, long first_row)
std::vector< T > createColumnData(const Euclid::Table::Table &table, size_t column_index)
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)
static std::string GenericScalarFormat(const Table &, size_t)
void populateNdArrayColumn(const Table &table, size_t column_index, CCfits::ExtHDU &table_hdu, long first_row)
T push_back(T... args)
T reserve(T... args)
T size(T... args)
T str(T... args)