CCfits  2.5
ColumnData.h
1 // Astrophysics Science Division,
2 // NASA/ Goddard Space Flight Center
3 // HEASARC
4 // http://heasarc.gsfc.nasa.gov
5 // e-mail: ccfits@legacy.gsfc.nasa.gov
6 //
7 // Original author: Ben Dorman
8 
9 #ifndef COLUMNDATA_H
10 #define COLUMNDATA_H 1
11 #include "CCfits.h"
12 
13 // vector
14 #include <vector>
15 // Column
16 #include "Column.h"
17 #ifdef _MSC_VER
18 #include "MSconfig.h"
19 #endif
20 
21 #include <complex>
22 #include <memory>
23 #include <iterator>
24 #include "FITSUtil.h"
25 using std::complex;
26 #include "FITS.h"
27 
28 
29 namespace CCfits {
30 
31 
32 
33  template <typename T>
34  class ColumnData : public Column //## Inherits: <unnamed>%385E51565EE8
35  {
36 
37  public:
38  ColumnData(const ColumnData< T > &right);
39  ColumnData (Table* p = 0);
40  ColumnData (int columnIndex, const string &columnName, ValueType type, const String &format, const String &unit, Table* p, int rpt = 1, long w = 1, const String &comment = "");
41  ~ColumnData();
42 
43  virtual ColumnData<T>* clone () const;
44  virtual void readData (long firstRow, long nelements, long firstElem = 1);
45  void setDataLimits (T* limits);
46  const T minLegalValue () const;
47  void minLegalValue (T value);
48  const T maxLegalValue () const;
49  void maxLegalValue (T value);
50  const T minDataValue () const;
51  void minDataValue (T value);
52  const T maxDataValue () const;
53  void maxDataValue (T value);
54  const std::vector<T>& data () const;
55  void setData (const std::vector<T>& value);
56  T data (int i);
57  void data (int i, T value);
58 
59  // Additional Public Declarations
60  friend class Column;
61  protected:
62  // Additional Protected Declarations
63 
64  private:
65  ColumnData< T > & operator=(const ColumnData< T > &right);
66 
67  void readColumnData (long firstRow, long nelements, T* nullValue = 0);
68  virtual bool compare (const Column &right) const;
69  virtual std::ostream& put (std::ostream& s) const;
70  void writeData (T* indata, long nRows = 1, long firstRow = 1, T* nullValue = 0);
71  void writeData (const std::vector<T>& indata, long firstRow = 1, T* nullValue = 0);
72  // Insert one or more blank rows into a FITS column.
73  virtual void insertRows (long first, long number = 1);
74  virtual void deleteRows (long first, long number = 1);
75 
76  // Additional Private Declarations
77 
78  private: //## implementation
79  // Data Members for Class Attributes
80  T m_minLegalValue;
81  T m_maxLegalValue;
82  T m_minDataValue;
83  T m_maxDataValue;
84 
85  // Data Members for Associations
86  std::vector<T> m_data;
87 
88  // Additional Implementation Declarations
89 
90  };
91 
92  // Parameterized Class CCfits::ColumnData
93 
94  template <typename T>
95  inline void ColumnData<T>::readData (long firstRow, long nelements, long firstElem)
96  {
97  readColumnData(firstRow,nelements,static_cast<T*>(0));
98  }
99 
100  template <typename T>
101  inline const T ColumnData<T>::minLegalValue () const
102  {
103  return m_minLegalValue;
104  }
105 
106  template <typename T>
107  inline void ColumnData<T>::minLegalValue (T value)
108  {
109  m_minLegalValue = value;
110  }
111 
112  template <typename T>
113  inline const T ColumnData<T>::maxLegalValue () const
114  {
115  return m_maxLegalValue;
116  }
117 
118  template <typename T>
119  inline void ColumnData<T>::maxLegalValue (T value)
120  {
121  m_maxLegalValue = value;
122  }
123 
124  template <typename T>
125  inline const T ColumnData<T>::minDataValue () const
126  {
127  return m_minDataValue;
128  }
129 
130  template <typename T>
131  inline void ColumnData<T>::minDataValue (T value)
132  {
133  m_minDataValue = value;
134  }
135 
136  template <typename T>
137  inline const T ColumnData<T>::maxDataValue () const
138  {
139  return m_maxDataValue;
140  }
141 
142  template <typename T>
143  inline void ColumnData<T>::maxDataValue (T value)
144  {
145  m_maxDataValue = value;
146  }
147 
148  template <typename T>
149  inline const std::vector<T>& ColumnData<T>::data () const
150  {
151  return m_data;
152  }
153 
154  template <typename T>
155  inline void ColumnData<T>::setData (const std::vector<T>& value)
156  {
157  m_data = value;
158  }
159 
160  template <typename T>
161  inline T ColumnData<T>::data (int i)
162  {
163  // return data stored in the ith row, which is in the i-1 th location in the array.
164  return m_data[i - 1];
165  }
166 
167  template <typename T>
168  inline void ColumnData<T>::data (int i, T value)
169  {
170  // assign data to i-1 th location in the array, representing the ith row.
171  m_data[i - 1] = value;
172  }
173 
174  // Parameterized Class CCfits::ColumnData
175 
176  template <typename T>
177  ColumnData<T>::ColumnData(const ColumnData<T> &right)
178  :Column(right),
179  m_minLegalValue(right.m_minLegalValue),
180  m_maxLegalValue(right.m_maxLegalValue),
181  m_minDataValue(right.m_minDataValue),
182  m_maxDataValue(right.m_maxDataValue),
183  m_data(right.m_data)
184  {
185  }
186 
187  template <typename T>
188  ColumnData<T>::ColumnData (Table* p)
189  : Column(p),
190  m_minLegalValue(),
191  m_maxLegalValue(),
192  m_minDataValue(),
193  m_maxDataValue(),
194  m_data()
195  {
196  }
197 
198  template <typename T>
199  ColumnData<T>::ColumnData (int columnIndex, const string &columnName, ValueType type, const String &format, const String &unit, Table* p, int rpt, long w, const String &comment)
200  : Column(columnIndex,columnName,type,format,unit,p,rpt,w,comment),
201  m_minLegalValue(),
202  m_maxLegalValue(),
203  m_minDataValue(),
204  m_maxDataValue(),
205  m_data()
206  {
207  }
208 
209 
210  template <typename T>
211  ColumnData<T>::~ColumnData()
212  {
213  }
214 
215 
216  template <typename T>
217  void ColumnData<T>::readColumnData (long firstRow, long nelements, T* nullValue)
218  {
219  if ( rows() < nelements )
220  {
221  std::cerr << "CCfits: More data requested than contained in table. ";
222  std::cerr << "Extracting complete column.\n";
223  nelements = rows();
224  }
225 
226  int status(0);
227  int anynul(0);
228 
229  FITSUtil::auto_array_ptr<T> array(new T[nelements]);
230 
231  makeHDUCurrent();
232 
233  if ( fits_read_col(fitsPointer(),type(), index(), firstRow, 1,
234  nelements, nullValue, array.get(), &anynul, &status) ) throw FitsError(status);
235 
236 
237  if (m_data.size() != static_cast<size_t>( rows() ) ) m_data.resize(rows());
238 
239  std::copy(&array[0],&array[nelements],m_data.begin()+firstRow-1);
240  if (nelements == rows()) isRead(true);
241  }
242 
243  template <typename T>
244  bool ColumnData<T>::compare (const Column &right) const
245  {
246  if ( !Column::compare(right) ) return false;
247  const ColumnData<T>& that = static_cast<const ColumnData<T>&>(right);
248  unsigned int n = m_data.size();
249  if ( that.m_data.size() != n ) return false;
250  for (unsigned int i = 0; i < n ; i++)
251  {
252  if (m_data[i] != that.m_data[i]) return false;
253  }
254  return true;
255  }
256 
257  template <typename T>
258  ColumnData<T>* ColumnData<T>::clone () const
259  {
260  return new ColumnData<T>(*this);
261  }
262 
263  template <typename T>
264  std::ostream& ColumnData<T>::put (std::ostream& s) const
265  {
266  Column::put(s);
267  if (FITS::verboseMode() && type() != Tstring)
268  {
269  s << " Column Legal limits: ( " << m_minLegalValue << "," << m_maxLegalValue << " )\n"
270  << " Column Data limits: ( " << m_minDataValue << "," << m_maxDataValue << " )\n";
271  }
272  if (!m_data.empty())
273  {
274  std::ostream_iterator<T> output(s,"\n");
275  // output each row on a separate line.
276  // user can supply manipulators to stream for formatting.
277  std::copy(m_data.begin(),m_data.end(),output);
278  }
279 
280  return s;
281  }
282 
283  template <typename T>
284  void ColumnData<T>::writeData (T* indata, long nRows, long firstRow, T* nullValue)
285  {
286 
287  // set columnData's data member to equal what's written to file.
288  // indata has size nRows: elements firstRow to firstRow + nRows - 1 will be written.
289  // if this exceeds the current rowlength of the HDU, update the return value for
290  // rows() in the parent after the fitsio call.
291  int status(0);
292  long elementsToWrite(nRows + firstRow -1);
293  // get a copy for restorative action.
294  std::vector<T> __tmp(m_data);
295 
296 
297  if (elementsToWrite > static_cast<long>(m_data.size()))
298  {
299 
300  m_data.resize(elementsToWrite,T());
301  }
302 
303  std::copy(&indata[0],&indata[nRows],m_data.begin()+firstRow-1);
304 
305  // if successful, write to disk.
306 
307  try
308  {
309  if (nullValue)
310  {
311  if (fits_write_colnull(fitsPointer(), type(), index(), firstRow, 1, nRows,
312  indata, nullValue, &status) != 0) throw FitsError(status);
313  }
314  else
315  {
316  if (fits_write_col(fitsPointer(), type(), index(), firstRow, 1, nRows,
317  indata, &status) != 0) throw FitsError(status);
318  }
319 
320  // tell the Table that the number of rows has changed
321  parent()->updateRows();
322  }
323  catch (FitsError) // the only thing that can throw here.
324  {
325  // reset to original content and rethrow the exception.
326  m_data = __tmp;
327  if (status == NO_NULL) throw NoNullValue(name());
328  else throw;
329  }
330  }
331 
332  template <typename T>
333  void ColumnData<T>::writeData (const std::vector<T>& indata, long firstRow, T* nullValue)
334  {
335  FITSUtil::CVarray<T> convert;
336  FITSUtil::auto_array_ptr<T> pcolData (convert(indata));
337  T* columnData = pcolData.get();
338  writeData(columnData,indata.size(),firstRow,nullValue);
339  }
340 
341  template <typename T>
342  void ColumnData<T>::insertRows (long first, long number)
343  {
344  FITSUtil::FitsNullValue<T> blank;
345  typename std::vector<T>::iterator in;
346  if (first !=0)
347  {
348  in = m_data.begin()+first;
349  }
350  else
351  {
352  in = m_data.begin();
353  }
354 
355  // non-throwing operations.
356  m_data.insert(in,number,blank());
357  }
358 
359  template <typename T>
360  void ColumnData<T>::deleteRows (long first, long number)
361  {
362  m_data.erase(m_data.begin()+first-1,m_data.begin()+first-1+number);
363  }
364 
365  template <typename T>
366  void ColumnData<T>::setDataLimits (T* limits)
367  {
368  m_minLegalValue = limits[0];
369  m_maxLegalValue = limits[1];
370  m_minDataValue = std::max(limits[2],limits[0]);
371  m_maxDataValue = std::min(limits[3],limits[1]);
372  }
373 
374  // Additional Declarations
375 
376  // all functions that operate on strings or complex data that call cfitsio
377  // need to be specialized.
378 
379 #if SPEC_TEMPLATE_IMP_DEFECT || SPEC_TEMPLATE_DECL_DEFECT
380 template <>
381 inline void ColumnData<complex<float> >::setDataLimits (complex<float>* limits)
382  {
383  m_minLegalValue = limits[0];
384  m_maxLegalValue = limits[1];
385  m_minDataValue = limits[2];
386  m_maxDataValue = limits[3];
387  }
388 #else
389 template <>
390  void ColumnData<complex<float> >::setDataLimits (complex<float>* limits);
391 #endif
392 
393 #if SPEC_TEMPLATE_IMP_DEFECT || SPEC_TEMPLATE_DECL_DEFECT
394 template <>
395 inline void ColumnData<complex<double> >::setDataLimits (complex<double>* limits)
396  {
397  m_minLegalValue = limits[0];
398  m_maxLegalValue = limits[1];
399  m_minDataValue = limits[2];
400  m_maxDataValue = limits[3];
401  }
402 #else
403  template <>
404  void ColumnData<complex<double> >::setDataLimits (complex<double>* limits);
405 #endif
406 
407 
408 #if SPEC_TEMPLATE_IMP_DEFECT || SPEC_TEMPLATE_DECL_DEFECT
409  template <>
410  inline void ColumnData<string>::readColumnData (long firstRow,
411  long nelements,
412  string* nullValue)
413  {
414  // nelements = nrows to read.
415  if (nelements < 1)
416  throw Column::InvalidNumberOfRows((int)nelements);
417  if (firstRow < 1 || (firstRow+nelements-1)>rows())
418  throw Column::InvalidRowNumber(name());
419 
420  int status = 0;
421  int anynul = 0;
422 
423  char** array = new char*[nelements];
424  // Initialize pointers to NULL so we can safely delete
425  // during error handling, even if they haven't been allocated.
426  for (long i=0; i<nelements; ++i)
427  array[i]=static_cast<char*>(0);
428  bool isError = false;
429 
430  // Strings are unusual. The variable length case is still
431  // handled by a ColumnData class, not a ColumnVectorData.
432  char* nulval = 0;
433  if (nullValue)
434  {
435  nulval = const_cast<char*>(nullValue->c_str());
436  }
437  else
438  {
439  nulval = new char;
440  *nulval = '\0';
441  }
442  makeHDUCurrent();
443  if (varLength())
444  {
445  long* strLengths = new long[nelements];
446  long* offsets = new long[nelements];
447  if (fits_read_descripts(fitsPointer(), index(), firstRow,
448  nelements, strLengths, offsets, &status))
449  {
450  isError = true;
451  }
452  else
453  {
454  // For variable length cols, must read 1 and only 1 row
455  // at a time into array.
456  for (long j=0; j<nelements; ++j)
457  {
458  array[j] = new char[strLengths[j] + 1];
459  }
460 
461  const long lastRow = firstRow+nelements-1;
462  for (long iRow=firstRow; !isError && iRow<=lastRow; ++iRow)
463  {
464  if (fits_read_col_str(fitsPointer(),index(), iRow, 1, 1,
465  nulval, &array[iRow-firstRow], &anynul,&status) )
466  isError=true;
467  }
468  }
469  delete [] strLengths;
470  delete [] offsets;
471  }
472  else
473  {
474  // Fixed length strings, length is stored in Column's m_width.
475  for (long j=0; j<nelements; ++j)
476  {
477  array[j] = new char[width() + 1];
478  }
479  if (fits_read_col_str(fitsPointer(),index(), firstRow,1,nelements,
480  nulval,array, &anynul,&status))
481  isError=true;
482  }
483 
484  if (isError)
485  {
486  // It's OK to do this even if error occurred before
487  // array rows were allocated. In that case their pointers
488  // were set to NULL.
489  for (long j = 0; j < nelements; ++j)
490  {
491  delete [] array[j];
492  }
493  delete [] array;
494  delete nulval;
495  throw FitsError(status);
496  }
497 
498  if (m_data.size() != static_cast<size_t>(rows()))
499  setData(std::vector<String>(rows(),String(nulval)));
500 
501  for (long j=0; j<nelements; ++j)
502  {
503  m_data[j - 1 + firstRow] = String(array[j]);
504  }
505 
506  for (long j=0; j<nelements; j++)
507  {
508  delete [] array[j];
509  }
510  delete [] array;
511  delete nulval;
512  if (nelements == rows()) isRead(true);
513 
514  }
515 #else
516  template <>
517 void ColumnData<string>::readColumnData (long firstRow, long nelements, string* nullValue);
518 #endif
519 
520 
521 #if SPEC_TEMPLATE_IMP_DEFECT || SPEC_TEMPLATE_DECL_DEFECT
522  template <>
523  inline void ColumnData<complex<float> >::readColumnData (long firstRow,
524  long nelements,
525  complex<float>* nullValue)
526  {
527  // specialization for ColumnData<string>
528  int status(0);
529  int anynul(0);
530  FITSUtil::auto_array_ptr<float> pArray(new float[nelements*2]);
531  float* array = pArray.get();
532  float nulval(0);
533  makeHDUCurrent();
534 
535 
536  if (fits_read_col_cmp(fitsPointer(),index(), firstRow,1,nelements,
537  nulval,array, &anynul,&status) ) throw FitsError(status);
538 
539 
540  if (m_data.size() != rows()) m_data.resize(rows());
541 
542  // the 'j -1 ' converts to zero based indexing.
543 
544  for (int j = 0; j < nelements; ++j)
545  {
546 
547  m_data[j - 1 + firstRow] = std::complex<float>(array[2*j],array[2*j+1]);
548  }
549  if (nelements == rows()) isRead(true);
550 
551  }
552 #else
553 template <>
554 void ColumnData<complex<float> >::readColumnData (long firstRow, long nelements,complex<float>* nullValue );
555 #endif
556 
557 #if SPEC_TEMPLATE_IMP_DEFECT || SPEC_TEMPLATE_DECL_DEFECT
558  template <>
559  inline void ColumnData<complex<double> >::readColumnData (long firstRow,
560  long nelements,
561  complex<double>* nullValue)
562  {
563  // specialization for ColumnData<complex<double> >
564  int status(0);
565  int anynul(0);
566  FITSUtil::auto_array_ptr<double> pArray(new double[nelements*2]);
567  double* array = pArray.get();
568  double nulval(0);
569  makeHDUCurrent();
570 
571 
572  if (fits_read_col_dblcmp(fitsPointer(), index(), firstRow,1,nelements,
573  nulval,array, &anynul,&status) ) throw FitsError(status);
574 
575 
576 
577 
578  if (m_data.size() != rows()) setData(std::vector<complex<double> >(rows(),nulval));
579 
580  // the 'j -1 ' converts to zero based indexing.
581 
582  for (int j = 0; j < nelements; j++)
583  {
584 
585  m_data[j - 1 + firstRow] = std::complex<double>(array[2*j],array[2*j+1]);
586  }
587  if (nelements == rows()) isRead(true);
588 
589  }
590 #else
591 template <>
592 void ColumnData<complex<double> >::readColumnData (long firstRow, long nelements,complex<double>* nullValue);
593 #endif
594 
595 #if SPEC_TEMPLATE_DECL_DEFECT
596  template <>
597  inline void ColumnData<string>::writeData (const std::vector<string>& indata,
598  long firstRow, string* nullValue)
599  {
600  int status=0;
601  char** columnData=FITSUtil::CharArray(indata);
602 
603  if ( fits_write_colnull(fitsPointer(), TSTRING, index(), firstRow, 1, indata.size(),
604  columnData, 0, &status) != 0 )
605  throw FitsError(status);
606  unsigned long elementsToWrite (indata.size() + firstRow - 1);
607  std::vector<string> __tmp(m_data);
608  if (m_data.size() < elementsToWrite)
609  {
610  m_data.resize(elementsToWrite,"");
611  std::copy(__tmp.begin(),__tmp.end(),m_data.begin());
612  }
613  std::copy(indata.begin(),indata.end(),m_data.begin()+firstRow-1);
614 
615 
616  for (size_t i = 0; i < indata.size(); ++i)
617  {
618  delete [] columnData[i];
619  }
620  delete [] columnData;
621  }
622 #else
623 template <>
624 void ColumnData<string>::writeData (const std::vector<string>& inData, long firstRow, string* nullValue);
625 #endif
626 
627 #ifdef SPEC_TEMPLATE_DECL_DEFECT
628  template <>
629  inline void ColumnData<complex<float> >::writeData (const std::vector<complex<float> >& inData,
630  long firstRow,
631  complex<float>* nullValue)
632  {
633  int status(0);
634  int nRows (inData.size());
635  FITSUtil::auto_array_ptr<float> pData(new float[nRows*2]);
636  float* Data = pData.get();
637  std::vector<complex<float> > __tmp(m_data);
638  for (int j = 0; j < nRows; ++j)
639  {
640  Data[ 2*j] = inData[j].real();
641  Data[ 2*j + 1] = inData[j].imag();
642  }
643 
644  try
645  {
646 
647  if (fits_write_col_cmp(fitsPointer(), index(), firstRow, 1,
648  nRows,Data, &status) != 0) throw FitsError(status);
649  long elementsToWrite(nRows + firstRow -1);
650  if (elementsToWrite > static_cast<long>(m_data.size()))
651  {
652 
653  m_data.resize(elementsToWrite);
654  }
655 
656  std::copy(inData.begin(),inData.end(),m_data.begin()+firstRow-1);
657 
658  // tell the Table that the number of rows has changed
659  parent()->updateRows();
660  }
661  catch (FitsError) // the only thing that can throw here.
662  {
663  // reset to original content and rethrow the exception.
664  m_data.resize(__tmp.size());
665  m_data = __tmp;
666  }
667 
668  }
669 
670 #else
671 template <>
672 void ColumnData<complex<float> >::writeData (const std::vector<complex<float> >& inData, long firstRow,
673  complex<float>* nullValue);
674 #endif
675 
676 #ifdef SPEC_TEMPLATE_DECL_DEFECT
677  template <>
678  inline void ColumnData<complex<double> >::writeData (const std::vector<complex<double> >& inData,
679  long firstRow,
680  complex<double>* nullValue)
681  {
682  int status(0);
683  int nRows (inData.size());
684  FITSUtil::auto_array_ptr<double> pData(new double[nRows*2]);
685  double* Data = pData.get();
686  std::vector<complex<double> > __tmp(m_data);
687  for (int j = 0; j < nRows; ++j)
688  {
689  pData[ 2*j] = inData[j].real();
690  pData[ 2*j + 1] = inData[j].imag();
691  }
692 
693  try
694  {
695 
696  if (fits_write_col_dblcmp(fitsPointer(), index(), firstRow, 1,
697  nRows,Data, &status) != 0) throw FitsError(status);
698  long elementsToWrite(nRows + firstRow -1);
699  if (elementsToWrite > static_cast<long>(m_data.size()))
700  {
701 
702  m_data.resize(elementsToWrite);
703  }
704 
705  std::copy(inData.begin(),inData.end(),m_data.begin()+firstRow-1);
706 
707  // tell the Table that the number of rows has changed
708  parent()->updateRows();
709  }
710  catch (FitsError) // the only thing that can throw here.
711  {
712  // reset to original content and rethrow the exception.
713  m_data.resize(__tmp.size());
714  m_data = __tmp;
715  }
716 
717  }
718 
719 #else
720 template <>
721 void ColumnData<complex<double> >::writeData (const std::vector<complex<double> >& inData, long firstRow,
722  complex<double>* nullValue);
723 
724 #endif
725 } // namespace CCfits
726 
727 
728 #endif
fitsfile * fitsPointer()
fits pointer corresponding to fits file containing column data.
Definition: Column.cxx:264
int index() const
get the Column index (the n in TTYPEn etc).
Definition: Column.h:1357
const String & format() const
return TFORMn keyword
Definition: Column.h:1519
Table * parent() const
return a pointer to the Table which owns this Column
Definition: Column.cxx:312
int rows() const
return the number of rows in the table.
Definition: Column.cxx:275
void updateRows()
update the number of rows in the table
Definition: Table.cxx:301
void makeHDUCurrent()
make HDU containing this the current HDU of the fits file.
Definition: Column.cxx:270
virtual std::ostream & put(std::ostream &s) const
internal implementation of << operator.
Definition: Column.cxx:302
static bool verboseMode()
return verbose setting for library
Definition: FITS.h:891
const String & comment() const
retrieve comment for Column
Definition: Column.h:1514
long width() const
return column data width
Definition: Column.h:1377
bool isRead() const
flag set to true if the entire column data has been read from disk
Definition: Column.h:1367
Namespace enclosing all CCfits classes and globals definitions.
Definition: AsciiTable.cxx:26
ValueType
CCfits value types and their CFITSIO equivalents (in caps)
Definition: CCfits.h:81
Column(const Column &right)
copy constructor, used in copying Columns to standard library containers.
Definition: Column.cxx:171
const String & unit() const
get units of data in Column (TUNITn keyword)
Definition: Column.h:1524
const String & name() const
return name of Column (TTYPEn keyword)
Definition: Column.h:1529
ValueType type() const
returns the data type of the column
Definition: Column.h:1434
bool varLength() const
boolean, set to true if Column has variable length vector rows.
Definition: Column.h:1392