Alexandria 2.25.0
SDC-CH common library for the Euclid project
FitsReaderHelper.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 "FitsReaderHelper.h"
27#include <CCfits/CCfits>
28#include <boost/lexical_cast.hpp>
29#include <boost/tokenizer.hpp>
30
31namespace Euclid {
32namespace Table {
33
35
36std::vector<std::string> autoDetectColumnNames(const CCfits::Table& table_hdu) {
38 for (int i = 1; i <= table_hdu.numCols(); ++i) {
39 std::string name = table_hdu.column(i).name();
40 if (name.empty()) {
41 name = "Col" + std::to_string(i);
42 }
43 names.push_back(std::move(name));
44 }
45 return names;
46}
47
49 if (format[0] == 'A') {
50 std::size_t size = std::stoi(format.substr(1));
51 return {typeid(std::string), size};
52 } else if (format[0] == 'I') {
53 return {typeid(int64_t), 0};
54 } else if (format[0] == 'F') {
55 return {typeid(double), 0};
56 } else if (format[0] == 'E') {
57 return {typeid(double), 0};
58 } else if (format[0] == 'D') {
59 return {typeid(double), 0};
60 }
61 throw Elements::Exception() << "FITS ASCII table format " << format << " is not "
62 << "yet supported";
63}
64
67 {'J', typeid(NdArray<int32_t>)}, {'B', typeid(NdArray<int32_t>)}, {'I', typeid(NdArray<int32_t>)},
68 {'K', typeid(NdArray<int64_t>)}, {'E', typeid(NdArray<float>)}, {'D', typeid(NdArray<double>)},
69 },
70 ScalarTypeMap{{'L', typeid(bool)}, {'J', typeid(int32_t)}, {'B', typeid(int32_t)}, {'I', typeid(int32_t)},
71 {'K', typeid(int64_t)}, {'E', typeid(float)}, {'D', typeid(double)}},
73 {'J', typeid(std::vector<int32_t>)}, {'K', typeid(std::vector<int64_t>)},
74 {'E', typeid(std::vector<float>)}, {'D', typeid(std::vector<double>)},
75 {'A', typeid(std::string)}};
76
78 const std::vector<size_t>& shape) {
79 // Get the size out of the format string
80 char ft = format.front();
81 int size = 1;
82 if (std::isdigit(format.front())) {
83 size = std::stoi(format.substr(0, format.size() - 1));
84 ft = format.back();
85 }
86
87 // If shape is set *and* it has more than one dimension, it is an NdArray
88 if (shape.size() > 1) {
89 auto i = std::find_if(NdTypeMap.begin(), NdTypeMap.end(),
90 [ft](const std::pair<char, std::type_index>& p) { return p.first == ft; });
91 if (i != NdTypeMap.end()) {
92 return {i->second, size};
93 }
94 }
95 // If the dimensionality is 1, it is a scalar
96 else if (size == 1) {
97 auto i = std::find_if(ScalarTypeMap.begin(), ScalarTypeMap.end(),
98 [ft](const std::pair<char, std::type_index>& p) { return p.first == ft; });
99 if (i != ScalarTypeMap.end()) {
100 return {i->second, size};
101 }
102 }
103 // Last, vectors
104 else {
105 auto i = std::find_if(VectorTypeMap.begin(), VectorTypeMap.end(),
106 [ft](const std::pair<char, std::type_index>& p) { return p.first == ft; });
107 if (i != VectorTypeMap.end()) {
108 return {i->second, size};
109 }
110 }
111 throw Elements::Exception() << "FITS binary table format " << format << " is not "
112 << "yet supported";
113}
114
116 std::vector<size_t> result{};
117 if (!tdim.empty() && tdim.front() == '(' && tdim.back() == ')') {
118 auto subtdim = tdim.substr(1, tdim.size() - 2);
119 boost::char_separator<char> sep{","};
120 boost::tokenizer<boost::char_separator<char>> tok{subtdim, sep};
121 for (auto& s : tok) {
122 result.push_back(boost::lexical_cast<size_t>(s));
123 }
124 // Note: the shape is in fortran order, so we need to reverse
125 std::reverse(result.begin(), result.end());
126 }
127 return result;
128}
129
132 for (int i = 1; i <= table_hdu.numCols(); i++) {
133 auto& column = table_hdu.column(i);
134
135 if (typeid(table_hdu) == typeid(CCfits::BinTable)) {
136 column.setDimen();
137 types.push_back(binaryFormatToType(column.format(), parseTDIM(column.dimen())));
138 } else {
139 types.push_back(asciiFormatToType(column.format()));
140 }
141 }
142 return types;
143}
144
145std::vector<std::string> autoDetectColumnUnits(const CCfits::Table& table_hdu) {
147 for (int i = 1; i <= table_hdu.numCols(); ++i) {
148 units.push_back(table_hdu.column(i).unit());
149 }
150 return units;
151}
152
154 std::vector<std::string> descriptions{};
155 for (int i = 1; i <= table_hdu.numCols(); ++i) {
156 std::string desc;
157 auto key = table_hdu.keyWord().find("TDESC" + std::to_string(i));
158 if (key != table_hdu.keyWord().end()) {
159 key->second->value(desc);
160 }
161 descriptions.push_back(desc);
162 }
163 return descriptions;
164}
165
166template <typename T>
167std::vector<Row::cell_type> convertScalarColumn(CCfits::Column& column, long first, long last) {
168 std::vector<T> data;
169 column.read(data, first, last);
171 for (auto value : data) {
172 result.push_back(value);
173 }
174 return result;
175}
176
177template <typename T>
178std::vector<Row::cell_type> convertVectorColumn(CCfits::Column& column, long first, long last) {
180 column.readArrays(data, first, last);
182 for (auto& valar : data) {
183 result.push_back(std::vector<T>(std::begin(valar), std::end(valar)));
184 }
185 return result;
186}
187
188template <typename T>
189std::vector<Row::cell_type> convertNdArrayColumn(CCfits::Column& column, long first, long last) {
191 column.readArrays(data, first, last);
192 std::vector<size_t> shape = parseTDIM(column.dimen());
193
195 for (auto& valar : data) {
196 result.push_back(NdArray<T>(shape, std::move(std::vector<T>(std::begin(valar), std::end(valar)))));
197 }
198 return result;
199}
200
202 return translateColumn(column, type, 1, column.rows());
203}
204
205std::vector<Row::cell_type> translateColumn(CCfits::Column& column, std::type_index type, long first, long last) {
206 if (type == typeid(bool)) {
207 return convertScalarColumn<bool>(column, first, last);
208 } else if (type == typeid(int32_t)) {
209 return convertScalarColumn<int32_t>(column, first, last);
210 } else if (type == typeid(int64_t)) {
211 return convertScalarColumn<int64_t>(column, first, last);
212 } else if (type == typeid(float)) {
213 return convertScalarColumn<float>(column, first, last);
214 } else if (type == typeid(double)) {
215 return convertScalarColumn<double>(column, first, last);
216 } else if (type == typeid(std::string)) {
217 return convertScalarColumn<std::string>(column, first, last);
218 } else if (type == typeid(std::vector<int32_t>)) {
219 return convertVectorColumn<int32_t>(column, first, last);
220 } else if (type == typeid(std::vector<int64_t>)) {
221 return convertVectorColumn<int64_t>(column, first, last);
222 } else if (type == typeid(std::vector<float>)) {
223 return convertVectorColumn<float>(column, first, last);
224 } else if (type == typeid(std::vector<double>)) {
225 return convertVectorColumn<double>(column, first, last);
226 } else if (type == typeid(NdArray<int32_t>)) {
227 return convertNdArrayColumn<int32_t>(column, first, last);
228 } else if (type == typeid(NdArray<int64_t>)) {
229 return convertNdArrayColumn<int64_t>(column, first, last);
230 } else if (type == typeid(NdArray<float>)) {
231 return convertNdArrayColumn<float>(column, first, last);
232 } else if (type == typeid(NdArray<double>)) {
233 return convertNdArrayColumn<double>(column, first, last);
234 }
235 throw Elements::Exception() << "Unsupported column type " << type.name();
236}
237
238} // namespace Table
239} // end of namespace Euclid
T back(T... args)
T begin(T... args)
NdArray(std::vector< size_t > shape_)
T empty(T... args)
T end(T... args)
T find_if(T... args)
T front(T... args)
T isdigit(T... args)
T move(T... args)
T name(T... args)
constexpr double s
std::vector< Row::cell_type > convertVectorColumn(CCfits::Column &column, long first, long last)
std::map< std::string, ColumnDescription > autoDetectColumnDescriptions(std::istream &in, const std::string &comment)
Reads the column descriptions of the given stream.
std::vector< Row::cell_type > convertNdArrayColumn(CCfits::Column &column, long first, long last)
std::pair< std::type_index, std::size_t > binaryFormatToType(const std::string &format, const std::vector< size_t > &shape)
std::pair< std::type_index, std::size_t > asciiFormatToType(const std::string &format)
const std::vector< std::pair< char, std::type_index > > ScalarTypeMap
const std::vector< std::pair< char, std::type_index > > NdTypeMap
std::vector< std::string > autoDetectColumnUnits(const CCfits::Table &table_hdu)
Reads the column units based on the TUNITn keyword.
std::vector< std::string > autoDetectColumnNames(std::istream &in, const std::string &comment, size_t columns_number)
Reads the column names of the given stream.
std::vector< Row::cell_type > translateColumn(CCfits::Column &column, std::type_index type)
Returns a vector representing the given FITS table column data, converted to the requested type.
std::vector< size_t > parseTDIM(const std::string &tdim)
const std::vector< std::pair< char, std::type_index > > VectorTypeMap
std::vector< std::pair< std::type_index, std::size_t > > autoDetectColumnTypes(const CCfits::Table &table_hdu)
Reads the column types of the given table HDU.
std::vector< Row::cell_type > convertScalarColumn(CCfits::Column &column, long first, long last)
T push_back(T... args)
T reverse(T... args)
T size(T... args)
T stoi(T... args)
T substr(T... args)
T to_string(T... args)