CCfits  2.6
ColumnVectorData.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 COLUMNVECTORDATA_H
10 #define COLUMNVECTORDATA_H 1
11 #ifdef _MSC_VER
12 #include "MSconfig.h"
13 #endif
14 #include "CCfits.h"
15 
16 // valarray
17 #include <valarray>
18 // vector
19 #include <vector>
20 // Column
21 #include "Column.h"
22 
23 #ifdef SSTREAM_DEFECT
24 #include <strstream>
25 #else
26 #include <sstream>
27 #endif
28 
29 #include <memory>
30 #include <numeric>
31 #include <algorithm>
32 
33 namespace CCfits {
34 
35  class Table;
36 
37 }
38 
39 #include "FITS.h"
40 #include "FITSUtil.h"
41 using std::complex;
42 
43 
44 namespace CCfits {
45 
46 
47 
48  template <typename T>
49  class ColumnVectorData : public Column //## Inherits: <unnamed>%38BAD1D4D370
50  {
51 
52  public:
53  ColumnVectorData(const ColumnVectorData< T > &right);
54  ColumnVectorData (Table* p = 0);
55  ColumnVectorData (int columnIndex, const string &columnName, ValueType type, const string &format, const string &unit, Table* p, int rpt = 1, long w = 1, const string &comment = "");
56  ~ColumnVectorData();
57 
58  virtual void readData (long firstrow, long nelements, long firstelem = 1);
59  virtual ColumnVectorData<T>* clone () const;
60  virtual void setDimen ();
61  void setDataLimits (T* limits);
62  const T minLegalValue () const;
63  void minLegalValue (T value);
64  const T maxLegalValue () const;
65  void maxLegalValue (T value);
66  const T minDataValue () const;
67  void minDataValue (T value);
68  const T maxDataValue () const;
69  void maxDataValue (T value);
70  const std::vector<std::valarray<T> >& data () const;
71  void setData (const std::vector<std::valarray<T> >& value);
72  const std::valarray<T>& data (int i) const;
73  void data (int i, const std::valarray<T>& value);
74 
75  // Additional Public Declarations
76  friend class Column;
77  protected:
78  // Additional Protected Declarations
79 
80  private:
81  ColumnVectorData< T > & operator=(const ColumnVectorData< T > &right);
82 
83  virtual bool compare (const Column &right) const;
84  void resizeDataObject (const std::vector<std::valarray<T> >& indata, size_t firstRow);
85  // Reads a specified number of column rows.
86  //
87  // There are no default arguments. The function
88  // Column::read(firstrow,firstelem,nelements)
89  // is designed for reading the whole column.
90  virtual void readColumnData (long first, long last, T* nullValue = 0);
91  virtual std::ostream& put (std::ostream& s) const;
92  void writeData (const std::valarray<T>& indata, long numRows, long firstRow = 1, T* nullValue = 0);
93  void writeData (const std::vector<std::valarray<T> >& indata, long firstRow = 1, T* nullValue = 0);
94  // Reads a specified number of column rows.
95  //
96  // There are no default arguments. The function
97  // Column::read(firstrow,firstelem,nelements)
98  // is designed for reading the whole column.
99  virtual void readRow (size_t row, T* nullValue = 0);
100  // Reads a variable row..
101  virtual void readVariableRow (size_t row, T* nullValue = 0);
102  void readColumnData (long firstrow, long nelements, long firstelem, T* nullValue = 0);
103  void writeData (const std::valarray<T>& indata, const std::vector<long>& vectorLengths, long firstRow = 1, T* nullValue = 0);
104  void writeFixedRow (const std::valarray<T>& data, long row, long firstElem = 1, T* nullValue = 0);
105  void writeFixedArray (T* data, long nElements, long nRows, long firstRow, T* nullValue = 0);
106  // Insert one or more blank rows into a FITS column.
107  virtual void insertRows (long first, long number = 1);
108  virtual void deleteRows (long first, long number = 1);
109  virtual size_t getStoredDataSize() const;
110  void doWrite (T* array, long row, long rowSize, long firstElem, T* nullValue);
111 
112  // Additional Private Declarations
113 
114  private: //## implementation
115  // Data Members for Class Attributes
116  T m_minLegalValue;
117  T m_maxLegalValue;
118  T m_minDataValue;
119  T m_maxDataValue;
120 
121  // Data Members for Associations
122  std::vector<std::valarray<T> > m_data;
123 
124  // Additional Implementation Declarations
125 
126  };
127 
128  // Parameterized Class CCfits::ColumnVectorData
129 
130  template <typename T>
131  inline void ColumnVectorData<T>::readData (long firstrow, long nelements, long firstelem)
132  {
133  readColumnData(firstrow,nelements,firstelem,static_cast<T*>(0));
134  }
135 
136  template <typename T>
137  inline const T ColumnVectorData<T>::minLegalValue () const
138  {
139  return m_minLegalValue;
140  }
141 
142  template <typename T>
143  inline void ColumnVectorData<T>::minLegalValue (T value)
144  {
145  m_minLegalValue = value;
146  }
147 
148  template <typename T>
149  inline const T ColumnVectorData<T>::maxLegalValue () const
150  {
151  return m_maxLegalValue;
152  }
153 
154  template <typename T>
155  inline void ColumnVectorData<T>::maxLegalValue (T value)
156  {
157  m_maxLegalValue = value;
158  }
159 
160  template <typename T>
161  inline const T ColumnVectorData<T>::minDataValue () const
162  {
163  return m_minDataValue;
164  }
165 
166  template <typename T>
167  inline void ColumnVectorData<T>::minDataValue (T value)
168  {
169  m_minDataValue = value;
170  }
171 
172  template <typename T>
173  inline const T ColumnVectorData<T>::maxDataValue () const
174  {
175  return m_maxDataValue;
176  }
177 
178  template <typename T>
179  inline void ColumnVectorData<T>::maxDataValue (T value)
180  {
181  m_maxDataValue = value;
182  }
183 
184  template <typename T>
185  inline const std::vector<std::valarray<T> >& ColumnVectorData<T>::data () const
186  {
187  return m_data;
188  }
189 
190  template <typename T>
191  inline void ColumnVectorData<T>::setData (const std::vector<std::valarray<T> >& value)
192  {
193  m_data = value;
194  }
195 
196  template <typename T>
197  inline const std::valarray<T>& ColumnVectorData<T>::data (int i) const
198  {
199  return m_data[i - 1];
200  }
201 
202  template <typename T>
203  inline void ColumnVectorData<T>::data (int i, const std::valarray<T>& value)
204  {
205  if (m_data[i-1].size() != value.size())
206  m_data[i-1].resize(value.size());
207  m_data[i - 1] = value;
208  }
209 
210  // Parameterized Class CCfits::ColumnVectorData
211 
212  template <typename T>
213  ColumnVectorData<T>::ColumnVectorData(const ColumnVectorData<T> &right)
214  :Column(right),
215  m_minLegalValue(right.m_minLegalValue),
216  m_maxLegalValue(right.m_maxLegalValue),
217  m_minDataValue(right.m_minDataValue),
218  m_maxDataValue(right.m_maxDataValue),
219  m_data(right.m_data)
220  {
221  }
222 
223  template <typename T>
224  ColumnVectorData<T>::ColumnVectorData (Table* p)
225  : Column(p),
226  m_minLegalValue(0),
227  m_maxLegalValue(0),
228  m_minDataValue(0),
229  m_maxDataValue(0),
230  m_data()
231  {
232  }
233 
234  template <typename T>
235  ColumnVectorData<T>::ColumnVectorData (int columnIndex, const string &columnName, ValueType type, const string &format, const string &unit, Table* p, int rpt, long w, const string &comment)
236  : Column(columnIndex,columnName,type,format,unit,p,rpt,w,comment),
237  m_minLegalValue(0),
238  m_maxLegalValue(0),
239  m_minDataValue(0),
240  m_maxDataValue(0),
241  m_data()
242  {
243  }
244 
245 
246  template <typename T>
247  ColumnVectorData<T>::~ColumnVectorData()
248  {
249  // valarray destructor should do all the work.
250  }
251 
252 
253  template <typename T>
254  bool ColumnVectorData<T>::compare (const Column &right) const
255  {
256  if ( !Column::compare(right) ) return false;
257  const ColumnVectorData<T>& that = static_cast<const ColumnVectorData<T>&>(right);
258  size_t n = m_data.size();
259  // m_data is of type valarray<T>.
260  if ( that.m_data.size() != n ) return false;
261  for (size_t i = 0; i < n ; i++)
262  {
263  const std::valarray<T>& thisValArray=m_data[i];
264  const std::valarray<T>& thatValArray=that.m_data[i];
265  size_t nn = thisValArray.size();
266  if (thatValArray.size() != nn ) return false;
267 
268  for (size_t j = 0; j < nn ; j++ )
269  {
270  if (thisValArray[j] != thatValArray[j])
271  return false;
272  }
273  }
274  return true;
275  }
276 
277  template <typename T>
278  ColumnVectorData<T>* ColumnVectorData<T>::clone () const
279  {
280  return new ColumnVectorData<T>(*this);
281  }
282 
283  template <typename T>
284  void ColumnVectorData<T>::resizeDataObject (const std::vector<std::valarray<T> >& indata, size_t firstRow)
285  {
286  // the rows() call is the value before updating.
287  // the updateRows() call at the end sets the call to return the
288  // value from the fits pointer - which is changed by writeFixedArray
289  // or writeFixedRow.
290 
291  const size_t lastInputRow(indata.size() + firstRow - 1);
292  const size_t newLastRow = std::max(lastInputRow,static_cast<size_t>(rows()));
293 
294  // if the write instruction increases the rows, we need to add
295  // rows to the data member and preserve its current contents.
296 
297  // rows() >= origNRows since it is the value for entire table,
298  // not just this column.
299  const size_t origNRows(m_data.size());
300  // This will always be an expansion. vector.resize() doesn't
301  // invalidate any data on an expansion.
302  if (newLastRow > origNRows) m_data.resize(newLastRow);
303 
304  if (varLength())
305  {
306  // The incoming data will determine each row size, thus
307  // no need to preserve any existing values in the row.
308  // Each value will eventually be overwritten.
309  for (size_t iRow = firstRow-1; iRow < lastInputRow; ++iRow)
310  {
311  std::valarray<T>& current = m_data[iRow];
312  const size_t newSize = indata[iRow - (firstRow-1)].size();
313  if (current.size() != newSize)
314  current.resize(newSize);
315  }
316  }
317  else
318  {
319  // All row sizes in m_data should ALWAYS be either repeat(),
320  // or 0 if they haven't been initialized. This is true regardless
321  // of the incoming data row size.
322 
323  // Perform LAZY initialization of m_data rows. Only
324  // expand a row valarray when it is first needed.
325  for (size_t iRow = firstRow-1; iRow < lastInputRow; ++iRow)
326  {
327  if (m_data[iRow].size() != repeat())
328  m_data[iRow].resize(repeat());
329  }
330  }
331  }
332 
333  template <typename T>
334  void ColumnVectorData<T>::setDimen ()
335  {
336  int status(0);
337  FITSUtil:: auto_array_ptr<char> dimValue (new char[FLEN_VALUE]);
338 
339 #ifdef SSTREAM_DEFECT
340  std::ostrstream key;
341 #else
342  std::ostringstream key;
343 #endif
344  key << "TDIM" << index();
345 
346 #ifdef SSTREAM_DEFECT
347  fits_read_key_str(fitsPointer(), key.str(), dimValue.get(),0,&status);
348 #else
349  fits_read_key_str(fitsPointer(),const_cast<char*>(key.str().c_str()),dimValue.get(),0,&status);
350 #endif
351 
352  if (status == 0)
353  {
354  dimen(String(dimValue.get()));
355  }
356  }
357 
358  template <typename T>
359  void ColumnVectorData<T>::readColumnData (long first, long last, T* nullValue)
360  {
361  makeHDUCurrent();
362 
363 
364  if ( rows() < last )
365  {
366  std::cerr << "CCfits: More data requested than contained in table. ";
367  std::cerr << "Extracting complete column.\n";
368  last = rows();
369  }
370 
371  long nelements = (last - first + 1)*repeat();
372 
373 
374  readColumnData(first,nelements,1,nullValue);
375  if (first <= 1 && last == rows()) isRead(true);
376  }
377 
378  template <typename T>
379  std::ostream& ColumnVectorData<T>::put (std::ostream& s) const
380  {
381  // output header information
382  Column::put(s);
383  if ( FITS::verboseMode() )
384  {
385  s << " Column Legal limits: ( " << m_minLegalValue << "," << m_maxLegalValue << " )\n"
386  << " Column Data limits: ( " << m_minDataValue << "," << m_maxDataValue << " )\n";
387  }
388  if (!m_data.empty())
389  {
390  for (size_t j = 0; j < m_data.size(); j++)
391  {
392  size_t n = m_data[j].size();
393  if ( n )
394  {
395  s << "Row " << j + 1 << " Vector Size " << n << '\n';
396  for (size_t k = 0; k < n - 1; k++)
397  {
398  s << m_data[j][k] << '\t';
399  }
400  s << m_data[j][n - 1] << '\n';
401  }
402  }
403  }
404 
405  return s;
406  }
407 
408  template <typename T>
409  void ColumnVectorData<T>::writeData (const std::valarray<T>& indata, long numRows, long firstRow, T* nullValue)
410  {
411  // This version of writeData is called by Column write functions which
412  // can only write the same number of elements to each row.
413  // For fixed width columns, this must be equal to the repeat value
414  // or an exception is thrown. For variable width, it only requires
415  // that indata.size()/numRows is an int.
416 
417  // won't do anything if < 0, and will give divide check if == 0.
418  if (numRows <= 0) throw InvalidNumberOfRows(numRows);
419 
420 #ifdef SSTREAM_DEFECT
421  std::ostrstream msgStr;
422 #else
423  std::ostringstream msgStr;
424 #endif
425  if (indata.size() % static_cast<size_t>(numRows))
426  {
427  msgStr << "To use this write function, input array size"
428  <<"\n must be exactly divisible by requested num rows: "
429  << numRows;
430  throw InsufficientElements(msgStr.str());
431  }
432  const size_t cellsize = indata.size()/static_cast<size_t>(numRows);
433 
434  if (!varLength() && cellsize != repeat() )
435  {
436  msgStr << "column: " << name()
437  << "\n input data size: " << indata.size()
438  << " required: " << numRows*repeat();
439  String msg(msgStr.str());
440  throw InsufficientElements(msg);
441  }
442 
443  std::vector<std::valarray<T> > internalFormat(numRows);
444 
445  // support writing equal row lengths to variable columns.
446 
447  for (long j = 0; j < numRows; ++j)
448  {
449  internalFormat[j].resize(cellsize);
450  internalFormat[j] = indata[std::slice(cellsize*j,cellsize,1)];
451  }
452 
453  // change the size of m_data based on the first row to be written
454  // and on the input data vector sizes.
455 
456  writeData(internalFormat,firstRow,nullValue);
457  }
458 
459  template <typename T>
460  void ColumnVectorData<T>::writeData (const std::vector<std::valarray<T> >& indata, long firstRow, T* nullValue)
461  {
462  // This is called directly by Column's writeArrays functions, and indirectly
463  // by both categories of write functions, ie. those which allow differing
464  // lengths per row and those that don't.
465  const size_t nInputRows(indata.size());
466  using std::valarray;
467 
468  resizeDataObject(indata,firstRow);
469  // After the above call, can assume all m_data arrays to be written to
470  // have been properly resized whether we're dealing with fixed or
471  // variable length.
472 
473  if (varLength())
474  {
475  // firstRow is 1-based, but all these internal row variables
476  // will be 0-based.
477  const size_t endRow = nInputRows + firstRow-1;
478  for (size_t iRow = firstRow-1; iRow < endRow; ++iRow)
479  {
480  m_data[iRow] = indata[iRow - (firstRow-1)];
481  // doWrite wants 1-based rows.
482  doWrite(&m_data[iRow][0], iRow+1, m_data[iRow].size(), 1, nullValue);
483  }
484  parent()->updateRows();
485  }
486  else
487  {
488  // Check for simplest case of all valarrays of size repeat().
489  // If any are greater, throw an error.
490  const size_t colRepeat = repeat();
491  bool allEqualRepeat = true;
492  for (size_t i=0; i<nInputRows; ++i)
493  {
494  const size_t sz = indata[i].size();
495  if (sz > colRepeat)
496  {
497 #ifdef SSTREAM_DEFECT
498  std::ostrstream oss;
499 #else
500  std::ostringstream oss;
501 #endif
502  oss << " vector column length " << colRepeat
503  <<", input valarray length " << sz;
504  throw InvalidRowParameter(oss.str());
505  }
506  if (sz < colRepeat)
507  allEqualRepeat = false;
508  }
509 
510  if (allEqualRepeat)
511  {
512  // concatenate the valarrays and write.
513  const size_t nElements (colRepeat*nInputRows);
514  FITSUtil::CVAarray<T> convert;
515  FITSUtil::auto_array_ptr<T> pArray(convert(indata));
516  T* array = pArray.get();
517 
518  // if T is complex, then CVAarray returns a
519  // C-array of complex objects. But FITS requires an array of complex's
520  // value_type.
521 
522  // This writes to the file and also calls updateRows.
523  writeFixedArray(array,nElements,nInputRows,firstRow,nullValue);
524 
525  for (size_t j = 0; j < nInputRows ; ++j)
526  {
527  const valarray<T>& input = indata[j];
528  valarray<T>& current = m_data[j + firstRow - 1];
529  // current should be resized by resizeDataObject.
530  current = input;
531  }
532  }
533  else
534  {
535  // Some input arrays have fewer than colRepeat elements.
536  const size_t endRow = nInputRows + firstRow-1;
537  for (size_t iRow = firstRow-1; iRow<endRow; ++iRow)
538  {
539  // resizeDataObject should already have resized all
540  // corresponding m_data rows to repeat().
541  const valarray<T>& input = indata[iRow-(firstRow-1)];
542  writeFixedRow(input, iRow, 1, nullValue);
543  }
544  parent()->updateRows();
545  }
546 
547  } // end if !varLength
548  }
549 
550  template <typename T>
551  void ColumnVectorData<T>::readRow (size_t row, T* nullValue)
552  {
553  makeHDUCurrent();
554 
555 
556 
557  if ( row > static_cast<size_t>(rows()) )
558  {
559 #ifdef SSTREAM_DEFECT
560  std::ostrstream msg;
561 #else
562  std::ostringstream msg;
563 #endif
564  msg << " row requested: " << row << " row range: 1 - " << rows();
565 #ifdef SSTREAM_DEFECT
566  msg << std::ends;
567 #endif
568 
569  throw Column::InvalidRowNumber(msg.str());
570  }
571 
572  // this is really for documentation purposes. I expect the optimizer will
573  // remove this redundant definition .
574  bool variable(type() < 0);
575 
576 
577  long nelements(repeat());
578 
579  if (variable)
580  {
581  readVariableRow(row,nullValue);
582  }
583  else
584  {
585  readColumnData(row,nelements,1,nullValue);
586  }
587  }
588 
589  template <typename T>
590  void ColumnVectorData<T>::readVariableRow (size_t row, T* nullValue)
591  {
592  int status(0);
593  long offset(0);
594  long repeat(0);
595  if (fits_read_descript(fitsPointer(),index(),static_cast<long>(row),
596  &repeat,&offset,&status)) throw FitsError(status);
597  readColumnData(row,repeat,1,nullValue);
598  }
599 
600  template <typename T>
601  void ColumnVectorData<T>::readColumnData (long firstrow, long nelements, long firstelem, T* nullValue)
602  {
603  int status=0;
604 
605  FITSUtil::auto_array_ptr<T> pArray(new T[nelements]);
606  T* array = pArray.get();
607  int anynul(0);
608 
609 
610 
611  if (fits_read_col(fitsPointer(), abs(type()),index(), firstrow, firstelem,
612  nelements, nullValue, array, &anynul, &status) != 0)
613  throw FitsError(status);
614 
615  size_t countRead = 0;
616  const size_t ONE = 1;
617 
618  if (m_data.size() != static_cast<size_t>(rows())) m_data.resize(rows());
619  size_t vectorSize(0);
620  if (!varLength())
621  {
622 
623  vectorSize = std::max(repeat(),ONE); // safety check.
624 
625  }
626  else
627  {
628  // assume that the user specified the correct length for
629  // variable columns. This should be ok since readVariableColumns
630  // uses fits_read_descripts to return this information from the
631  // fits pointer, and this is passed as nelements here.
632  vectorSize = nelements;
633  }
634  size_t n = nelements;
635 
636  int i = firstrow;
637  int ii = i - 1;
638  while ( countRead < n)
639  {
640  std::valarray<T>& current = m_data[ii];
641  if (current.size() != vectorSize) current.resize(vectorSize);
642  int elementsInFirstRow = vectorSize-firstelem + 1;
643  bool lastRow = ( (nelements - countRead) < vectorSize);
644  if (lastRow)
645  {
646  int elementsInLastRow = nelements - countRead;
647  std::valarray<T> ttmp(array + vectorSize*(ii-firstrow) + elementsInFirstRow,
648  elementsInLastRow);
649  for (int kk = 0; kk < elementsInLastRow; kk++) current[kk] = ttmp[kk];
650  countRead += elementsInLastRow;
651 
652  }
653  // what to do with complete rows
654  else
655  {
656  if (firstelem == 1 || (firstelem > 1 && i > firstrow) )
657  {
658  std::valarray<T> ttmp(array + vectorSize*(ii - firstrow) +
659  elementsInFirstRow,vectorSize);
660  current = ttmp;
661  ii++;
662  i++;
663  countRead += vectorSize;
664  }
665  else
666  {
667  if (i == firstrow)
668  {
669  std::valarray<T> ttmp(array,elementsInFirstRow);
670  for (size_t kk = firstelem ; kk < vectorSize ; kk++)
671  current[kk] = ttmp[kk-firstelem];
672  countRead += elementsInFirstRow;
673  i++;
674  ii++;
675  }
676  }
677  }
678  }
679  }
680 
681  template <typename T>
682  void ColumnVectorData<T>::writeData (const std::valarray<T>& indata, const std::vector<long>& vectorLengths, long firstRow, T* nullValue)
683  {
684  // Called from Column write functions which allow differing lengths
685  // for each row.
686  using namespace std;
687  const size_t N(vectorLengths.size());
688  vector<long> sums(N);
689  // pre-calculate partial sums of vector lengths for use as array offsets.
690  partial_sum(vectorLengths.begin(),vectorLengths.end(),sums.begin());
691  // check that sufficient data have been supplied to carry out the entire operation.
692  if (indata.size() < static_cast<size_t>(sums[N-1]) )
693  {
694 #ifdef SSTREAM_DEFECT
695  ostrstream msgStr;
696 #else
697  ostringstream msgStr;
698 #endif
699  msgStr << " input data size: " << indata.size() << " vector length sum: " << sums[N-1];
700 #ifdef SSTREAM_DEFECT
701  msgStr << std::ends;
702 #endif
703 
704  String msg(msgStr.str());
705  throw InsufficientElements(msg);
706  }
707 
708  vector<valarray<T> > vvArray(N);
709  long& last = sums[0];
710  vvArray[0].resize(last);
711  for (long jj = 0; jj < last; ++jj) vvArray[0][jj] = indata[jj];
712 
713  for (size_t j = 1; j < N; ++j)
714  {
715  valarray<T>& __tmp = vvArray[j];
716  // these make the code much more readable
717  long& first = sums[j-1];
718  long& jlast = sums[j];
719  __tmp.resize(jlast - first);
720  for (long k = first; k < jlast; ++k)
721  {
722  __tmp[k - first] = indata[k];
723  }
724  }
725 
726  writeData(vvArray,firstRow,nullValue);
727  }
728 
729  template <typename T>
730  void ColumnVectorData<T>::writeFixedRow (const std::valarray<T>& data, long row, long firstElem, T* nullValue)
731  {
732 
733  // This is to be called only for FIXED length vector columns. It will
734  // throw if data.size()+firstElem goes beyond the repeat value.
735  // If data.size() is less than repeat, it leaves the remaining values
736  // undisturbed both in the file and in m_data storage.
737 
738 #ifdef SSTREAM_DEFECT
739  std::ostrstream msgStr;
740 #else
741  std::ostringstream msgStr;
742 #endif
743  if (varLength())
744  {
745  msgStr <<"Calling ColumnVectorData::writeFixedRow for a variable length column.\n";
746  throw FitsFatal(msgStr.str());
747  }
748 
749  std::valarray<T>& storedRow = m_data[row];
750  long inputSize = static_cast<long>(data.size());
751  long storedSize(storedRow.size());
752  if (storedSize != static_cast<long>(repeat()))
753  {
754  msgStr<<"stored array size vs. column width mismatch in ColumnVectorData::writeFixedRow.\n";
755  throw FitsFatal(msgStr.str());
756  }
757 
758  if (inputSize + firstElem - 1 > storedSize)
759  {
760  msgStr << " requested write " << firstElem << " to "
761  << firstElem + inputSize - 1 << " exceeds vector length " << repeat();
762  throw InvalidRowParameter(msgStr.str());
763  }
764 
765  // CANNOT give a strong exception safety guarantee because writing
766  // data changes the file. Any corrective action that could be taken
767  // [e.g. holding initial contents of the row and writing it back after
768  // an exception is thrown] could in principle throw the same exception
769  // we are trying to protect from.
770 
771  // routine does however give the weak guarantee (no resource leaks).
772 
773  // It's never a good thing to cast away a const, but doWrite calls the
774  // CFITSIO write functions which take a non-const pointer (though
775  // it shouldn't actually modify the array), and I'd rather not
776  // copy the entire valarray just to avoid this problem.
777  std::valarray<T>& lvData = const_cast<std::valarray<T>&>(data);
778  T* inPointer = &lvData[0];
779  doWrite(inPointer, row+1, inputSize, firstElem, nullValue);
780 
781  // Writing to disk was successful, now update FITS object and return.
782  const size_t offset = static_cast<size_t>(firstElem) - 1;
783  for (size_t iElem=0; iElem < static_cast<size_t>(inputSize); ++iElem)
784  {
785  // This doesn't require inPointer's non-constness. It's just
786  // used here to speed things up a bit.
787  storedRow[iElem + offset] = inPointer[iElem];
788  }
789  }
790 
791  template <typename T>
792  void ColumnVectorData<T>::writeFixedArray (T* data, long nElements, long nRows, long firstRow, T* nullValue)
793  {
794  int status(0);
795 
796  // check for sanity of inputs, then write to file.
797  // this function writes only complete rows to a table with
798  // fixed width rows.
799 
800 
801  if ( nElements < nRows*static_cast<long>(repeat()) )
802  {
803 #ifdef SSTREAM_DEFECT
804  std::ostrstream msgStr;
805 #else
806  std::ostringstream msgStr;
807 #endif
808  msgStr << " input array size: " << nElements << " required " << nRows*repeat();
809  String msg(msgStr.str());
810 
811  throw Column::InsufficientElements(msg);
812  }
813 
814  if (nullValue)
815  {
816  if (fits_write_colnull(fitsPointer(),abs(type()),index(),firstRow,
817  1,nElements,data,nullValue,&status)) throw FitsError(status);
818  }
819  else
820  {
821  if (fits_write_col(fitsPointer(),abs(type()),index(),firstRow,
822  1,nElements,data,&status)) throw FitsError(status);
823  }
824 
825  parent()->updateRows();
826  }
827 
828  template <typename T>
829  void ColumnVectorData<T>::insertRows (long first, long number)
830  {
831  if (first >= 0 && first <= static_cast<long>(m_data.size()))
832  {
833  typename std::vector<std::valarray<T> >::iterator in;
834  if (first !=0)
835  {
836  in = m_data.begin()+first;
837  }
838  else
839  {
840  in = m_data.begin();
841  }
842 
843  // non-throwing operations.
844  m_data.insert(in,number,std::valarray<T>(T(),0));
845  }
846  }
847 
848  template <typename T>
849  void ColumnVectorData<T>::deleteRows (long first, long number)
850  {
851  // Don't assume the calling routine (ie. Table's deleteRows)
852  // knows Column's current m_data size. m_data may still be
853  // size 0 if no read operations have been performed on Column.
854  // Therefore perform range checking before erasing.
855  const long curSize = static_cast<long>(m_data.size());
856  if (curSize>0 && first <= curSize)
857  {
858  const long last = std::min(curSize, first-1+number);
859  m_data.erase(m_data.begin()+first-1,m_data.begin()+last);
860  }
861 
862  }
863 
864  template <typename T>
865  size_t ColumnVectorData<T>::getStoredDataSize() const
866  {
867  return m_data.size();
868  }
869 
870  template <typename T>
871  void ColumnVectorData<T>::setDataLimits (T* limits)
872  {
873  m_minLegalValue = limits[0];
874  m_maxLegalValue = limits[1];
875  m_minDataValue = std::max(limits[2],limits[0]);
876  m_maxDataValue = std::min(limits[3],limits[1]);
877  }
878 
879  template <typename T>
880  void ColumnVectorData<T>::doWrite (T* array, long row, long rowSize, long firstElem, T* nullValue)
881  {
882  int status(0);
883  // internal functioning of write_colnull forbids its use for writing
884  // variable width columns. If a nullvalue argument was supplied it will
885  // be ignored.
886  if ( !varLength())
887  {
888  if (fits_write_colnull(fitsPointer(),type(),index(),row, firstElem, rowSize,
889  array, nullValue,&status)) throw FitsError(status);
890  }
891  else
892  {
893  if (fits_write_col(fitsPointer(),abs(type()),index(),row,firstElem,rowSize,
894  array,&status)) throw FitsError(status);
895 
896  }
897  }
898 
899  // Additional Declarations
900 
901  // all functions that operate on complex data that call cfitsio
902  // need to be specialized. The signature containing complex<T>* objects
903  // is unfortunate, perhaps, for this purpose, but the user will access
904  // rw operations through standard library containers.
905 
906 
907 
908 
909 
910 #if SPEC_TEMPLATE_IMP_DEFECT || SPEC_TEMPLATE_DECL_DEFECT
911 template <>
912 inline void ColumnVectorData<complex<float> >::setDataLimits (complex<float>* limits)
913  {
914  m_minLegalValue = limits[0];
915  m_maxLegalValue = limits[1];
916  m_minDataValue = limits[2];
917  m_maxDataValue = limits[3];
918  }
919 #else
920 template <>
921  void
922  ColumnVectorData<complex<float> >::setDataLimits (complex<float>* limits);
923 #endif
924 
925 #if SPEC_TEMPLATE_IMP_DEFECT || SPEC_TEMPLATE_DECL_DEFECT
926 template <>
927 inline void ColumnVectorData<complex<double> >::setDataLimits (complex<double>* limits)
928  {
929  m_minLegalValue = limits[0];
930  m_maxLegalValue = limits[1];
931  m_minDataValue = limits[2];
932  m_maxDataValue = limits[3];
933  }
934 #else
935  template <>
936  void
937  ColumnVectorData<complex<double> >::setDataLimits (complex<double>* limits);
938 #endif
939 
940 
941 #if SPEC_TEMPLATE_IMP_DEFECT || SPEC_TEMPLATE_DECL_DEFECT
942  template <>
943  inline void ColumnVectorData<std::complex<float> >::readColumnData(long firstRow,
944  long nelements, long firstElem, std::complex<float>* null )
945  {
946  int status=0;
947  float nulval (0);
948  FITSUtil::auto_array_ptr<float> pArray(new float[2*nelements]);
949  float* array = pArray.get();
950  int anynul(0);
951 
952  if (fits_read_col_cmp(fitsPointer(),index(),firstRow, firstElem,
953  nelements,nulval,array,&anynul,&status) ) throw FitsError(status);
954 
955  if (m_data.size() != static_cast<size_t>(rows())) m_data.resize(rows());
956 
957  std::valarray<std::complex<float> > readData(nelements);
958  for (long j = 0; j < nelements; ++j)
959  {
960  readData[j] = std::complex<float>(array[2*j],array[2*j+1]);
961  }
962  size_t countRead = 0;
963  const size_t ONE = 1;
964 
965  if (m_data.size() != static_cast<size_t>(rows())) m_data.resize(rows());
966  size_t vectorSize(0);
967  if (!varLength())
968  {
969  vectorSize = std::max(repeat(),ONE); // safety check.
970  }
971  else
972  {
973  // assume that the user specified the correct length for
974  // variable columns. This should be ok since readVariableColumns
975  // uses fits_read_descripts to return this information from the
976  // fits pointer, and this is passed as nelements here.
977  vectorSize = nelements;
978  }
979  size_t n = nelements;
980 
981  int i = firstRow;
982  int ii = i - 1;
983  while ( countRead < n)
984  {
985  std::valarray<complex<float> >& current = m_data[ii];
986  if (current.size() != vectorSize) current.resize(vectorSize,0.);
987  int elementsInFirstRow = vectorSize-firstElem + 1;
988  bool lastRow = ( (nelements - countRead) < vectorSize);
989  if (lastRow)
990  {
991  int elementsInLastRow = nelements - countRead;
992  std::copy(&readData[countRead],&readData[0]+nelements,&current[0]);
993  countRead += elementsInLastRow;
994  }
995  // what to do with complete rows. if firstElem == 1 the
996  else
997  {
998  if (firstElem == 1 || (firstElem > 1 && i > firstRow) )
999  {
1000  current = readData[std::slice(vectorSize*(ii-firstRow)+
1001  elementsInFirstRow,vectorSize,1)];
1002  ++ii;
1003  ++i;
1004  countRead += vectorSize;
1005  }
1006  else
1007  {
1008  if (i == firstRow)
1009  {
1010  std::copy(&readData[0],&readData[0]+elementsInFirstRow,
1011  &current[firstElem]);
1012  countRead += elementsInFirstRow;
1013  ++i;
1014  ++ii;
1015  }
1016  }
1017  }
1018  }
1019  }
1020 #else
1021 template <>
1022 void ColumnVectorData<complex<float> >::readColumnData(long firstRow,
1023  long nelements,
1024  long firstElem, complex<float>* null);
1025 #endif
1026 
1027 #if SPEC_TEMPLATE_IMP_DEFECT || SPEC_TEMPLATE_DECL_DEFECT
1028  template <>
1029  inline void ColumnVectorData<complex<double> >::readColumnData (long firstRow,
1030  long nelements,long firstElem,
1031  complex<double>* nullValue)
1032  {
1033 
1034  // duplicated for each complex type to work around imagined or
1035  // actual compiler deficiencies.
1036  int status=0;
1037  double nulval (0);
1038  FITSUtil::auto_array_ptr<double> pArray(new double[2*nelements]);
1039  double* array = pArray.get();
1040  int anynul(0);
1041 
1042  if (fits_read_col_dblcmp(fitsPointer(),index(),firstRow, firstElem,
1043  nelements,nulval,array,&anynul,&status) ) throw FitsError(status);
1044 
1045  if (m_data.size() != static_cast<size_t>(rows())) m_data.resize(rows());
1046 
1047  std::valarray<std::complex<double> > readData(nelements);
1048  for (long j = 0; j < nelements; ++j)
1049  {
1050  readData[j] = std::complex<double>(array[2*j],array[2*j+1]);
1051  }
1052  size_t countRead = 0;
1053  const size_t ONE = 1;
1054 
1055  if (m_data.size() != static_cast<size_t>(rows())) m_data.resize(rows());
1056  size_t vectorSize(0);
1057  if (!varLength())
1058  {
1059  vectorSize = std::max(repeat(),ONE); // safety check.
1060  }
1061  else
1062  {
1063  // assume that the user specified the correct length for
1064  // variable columns. This should be ok since readVariableColumns
1065  // uses fits_read_descripts to return this information from the
1066  // fits pointer, and this is passed as nelements here.
1067  vectorSize = nelements;
1068  }
1069  size_t n = nelements;
1070 
1071  int i = firstRow;
1072  int ii = i - 1;
1073  while ( countRead < n)
1074  {
1075  std::valarray<std::complex<double> >& current = m_data[ii];
1076  if (current.size() != vectorSize) current.resize(vectorSize,0.);
1077  int elementsInFirstRow = vectorSize-firstElem + 1;
1078  bool lastRow = ( (nelements - countRead) < vectorSize);
1079  if (lastRow)
1080  {
1081  int elementsInLastRow = nelements - countRead;
1082  std::copy(&readData[countRead],&readData[0]+nelements,&current[0]);
1083  countRead += elementsInLastRow;
1084  }
1085  // what to do with complete rows. if firstElem == 1 the
1086  else
1087  {
1088  if (firstElem == 1 || (firstElem > 1 && i > firstRow) )
1089  {
1090  current = readData[std::slice(vectorSize*(ii-firstRow)+
1091  elementsInFirstRow,vectorSize,1)];
1092  ++ii;
1093  ++i;
1094  countRead += vectorSize;
1095  }
1096  else
1097  {
1098  if (i == firstRow)
1099  {
1100  std::copy(&readData[0],&readData[0]+elementsInFirstRow,
1101  &current[firstElem]);
1102  countRead += elementsInFirstRow;
1103  ++i;
1104  ++ii;
1105  }
1106  }
1107  }
1108  }
1109  }
1110 #else
1111 template <>
1112 void ColumnVectorData<complex<double> >::readColumnData (long firstRow,
1113  long nelements,
1114  long firstElem, complex<double>* null);
1115 #endif
1116 
1117 #if SPEC_TEMPLATE_IMP_DEFECT || SPEC_TEMPLATE_DECL_DEFECT
1118  template <>
1119  inline void ColumnVectorData<complex<float> >::writeFixedArray
1120  (complex<float>* data, long nElements, long nRows, long firstRow,
1121  complex<float>* nullValue)
1122  {
1123 
1124  int status(0);
1125 
1126  // check for sanity of inputs, then write to file.
1127  // this function writes only complete rows to a table with
1128  // fixed width rows.
1129 
1130 
1131  if ( nElements < nRows*static_cast<long>(repeat()) )
1132  {
1133 #ifdef SSTREAM_DEFECT
1134  std::ostrstream msgStr;
1135 #else
1136  std::ostringstream msgStr;
1137 #endif
1138  msgStr << " input array size: " << nElements
1139  << " required " << nRows*repeat();
1140 #ifdef SSTREAM_DEFECT
1141  msgStr << std::ends;
1142 #endif
1143 
1144 
1145  String msg(msgStr.str());
1146 
1147  throw Column::InsufficientElements(msg);
1148  }
1149 
1150  FITSUtil::auto_array_ptr<float> realData(new float[2*nElements]);
1151 
1152  for (int j = 0; j < nElements; ++j)
1153  {
1154  realData[2*j] = data[j].real();
1155  realData[2*j+1] = data[j].imag();
1156  }
1157 
1158 
1159 
1160  if (fits_write_col_cmp(fitsPointer(),index(),firstRow,
1161  1,nElements,realData.get(),&status)) throw FitsError(status);
1162 
1163  parent()->updateRows();
1164  }
1165 #else
1166 template <>
1167 void ColumnVectorData<complex<float> >::writeFixedArray
1168  (complex<float>* data, long nElements, long nRows, long firstRow, std::complex<float>* null);
1169 #endif
1170 
1171 #if SPEC_TEMPLATE_IMP_DEFECT || SPEC_TEMPLATE_DECL_DEFECT
1172  template <>
1173  inline void ColumnVectorData<complex<double> >::writeFixedArray
1174  (complex<double>* data, long nElements, long nRows, long firstRow,
1175  complex<double>* nullValue)
1176  {
1177  int status(0);
1178 
1179  // check for sanity of inputs, then write to file.
1180  // this function writes only complete rows to a table with
1181  // fixed width rows.
1182 
1183 
1184  if ( nElements < nRows*static_cast<long>(repeat()) )
1185  {
1186 #ifdef SSTREAM_DEFECT
1187  std::ostrstream msgStr;
1188 #else
1189  std::ostringstream msgStr;
1190 #endif
1191  msgStr << " input array size: " << nElements
1192  << " required " << nRows*repeat();
1193 #ifdef SSTREAM_DEFECT
1194  msgStr << std::ends;
1195 #endif
1196 
1197  String msg(msgStr.str());
1198 
1199  throw Column::InsufficientElements(msg);
1200  }
1201 
1202  FITSUtil::auto_array_ptr<double> realData(new double[2*nElements]);
1203 
1204  for (int j = 0; j < nElements; ++j)
1205  {
1206  realData[2*j] = data[j].real();
1207  realData[2*j+1] = data[j].imag();
1208  }
1209 
1210 
1211 
1212  if (fits_write_col_dblcmp(fitsPointer(),index(),firstRow,
1213  1,nElements,realData.get(),&status)) throw FitsError(status);
1214 
1215  parent()->updateRows();
1216 
1217  }
1218 #else
1219 template <>
1220 void ColumnVectorData<complex<double> >::writeFixedArray
1221  (complex<double>* data, long nElements, long nRows, long firstRow,
1222  std::complex<double>* null);
1223 #endif
1224 
1225 #ifdef SPEC_TEMPLATE_DECL_DEFECT
1226  template <>
1227  inline void
1228  ColumnVectorData<std::complex<float> >::doWrite
1229  (std::complex<float>* data, long row, long rowSize, long firstElem, std::complex<float>* nullValue )
1230  {
1231  int status(0);
1232  FITSUtil::auto_array_ptr<float> carray( new float[2*rowSize]);
1233  for ( long j = 0 ; j < rowSize; ++ j)
1234  {
1235  carray[2*j] = data[j].real();
1236  carray[2*j + 1] = data[j].imag();
1237  }
1238  if (fits_write_col_cmp(fitsPointer(),index(),row,firstElem,rowSize,
1239  carray.get(),&status)) throw FitsError(status);
1240  }
1241 
1242 
1243  template <>
1244  inline void
1245  ColumnVectorData<std::complex<double> >::doWrite
1246  (std::complex<double>* data, long row, long rowSize, long firstElem, std::complex<double>* nullValue )
1247  {
1248  int status(0);
1249  FITSUtil::auto_array_ptr<double> carray( new double[2*rowSize]);
1250  for ( long j = 0 ; j < rowSize; ++ j)
1251  {
1252  carray[2*j] = data[j].real();
1253  carray[2*j + 1] = data[j].imag();
1254  }
1255  if (fits_write_col_dblcmp(fitsPointer(),index(),row,firstElem,rowSize,
1256  carray.get(),&status)) throw FitsError(status);
1257 
1258  }
1259 
1260 #else
1261 template<>
1262 void
1263 ColumnVectorData<complex<float> >::doWrite
1264  ( complex<float>* data, long row, long rowSize, long firstElem, complex<float>* nullValue);
1265 
1266 template<>
1267 void
1268 ColumnVectorData<complex<double> >::doWrite
1269  ( complex<double>* data, long row, long rowSize, long firstElem, complex<double>* nullValue );
1270 #endif
1271 } // namespace CCfits
1272 
1273 
1274 #endif
virtual std::ostream & put(std::ostream &s) const
internal implementation of &lt;&lt; operator.
Definition: Column.cxx:302
static bool verboseMode()
return verbose setting for library
Definition: FITS.h:862
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