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