00001
00002
00003
00004
00005
00006
00007
00008
00009 #ifndef COLUMNDATA_H
00010 #define COLUMNDATA_H 1
00011 #include "CCfits.h"
00012
00013
00014 #include <vector>
00015
00016 #include "Column.h"
00017 #ifdef _MSC_VER
00018 #include "MSconfig.h"
00019 #endif
00020
00021 #include <complex>
00022 #include <memory>
00023 #include <iterator>
00024 #include "FITSUtil.h"
00025 using std::complex;
00026 #include "FITS.h"
00027
00028
00029 namespace CCfits {
00030
00031
00032
00033 template <typename T>
00034 class ColumnData : public Column
00035 {
00036
00037 public:
00038 ColumnData(const ColumnData< T > &right);
00039 ColumnData (Table* p = 0);
00040 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 = "");
00041 ~ColumnData();
00042
00043 virtual ColumnData<T>* clone () const;
00044 virtual void readData (long firstRow, long nelements, long firstElem = 1);
00045 void setDataLimits (T* limits);
00046 const T minLegalValue () const;
00047 void minLegalValue (T value);
00048 const T maxLegalValue () const;
00049 void maxLegalValue (T value);
00050 const T minDataValue () const;
00051 void minDataValue (T value);
00052 const T maxDataValue () const;
00053 void maxDataValue (T value);
00054 const std::vector<T>& data () const;
00055 void setData (const std::vector<T>& value);
00056 T data (int i);
00057 void data (int i, T value);
00058
00059
00060 friend class Column;
00061 protected:
00062
00063
00064 private:
00065 ColumnData< T > & operator=(const ColumnData< T > &right);
00066
00067 void readColumnData (long firstRow, long nelements, T* nullValue = 0);
00068 virtual bool compare (const Column &right) const;
00069 virtual std::ostream& put (std::ostream& s) const;
00070 void writeData (T* indata, long nRows = 1, long firstRow = 1, T* nullValue = 0);
00071 void writeData (const std::vector<T>& indata, long firstRow = 1, T* nullValue = 0);
00072
00073 virtual void insertRows (long first, long number = 1);
00074 virtual void deleteRows (long first, long number = 1);
00075
00076
00077
00078 private:
00079
00080 T m_minLegalValue;
00081 T m_maxLegalValue;
00082 T m_minDataValue;
00083 T m_maxDataValue;
00084
00085
00086 std::vector<T> m_data;
00087
00088
00089
00090 };
00091
00092
00093
00094 template <typename T>
00095 inline void ColumnData<T>::readData (long firstRow, long nelements, long firstElem)
00096 {
00097 readColumnData(firstRow,nelements,static_cast<T*>(0));
00098 }
00099
00100 template <typename T>
00101 inline const T ColumnData<T>::minLegalValue () const
00102 {
00103 return m_minLegalValue;
00104 }
00105
00106 template <typename T>
00107 inline void ColumnData<T>::minLegalValue (T value)
00108 {
00109 m_minLegalValue = value;
00110 }
00111
00112 template <typename T>
00113 inline const T ColumnData<T>::maxLegalValue () const
00114 {
00115 return m_maxLegalValue;
00116 }
00117
00118 template <typename T>
00119 inline void ColumnData<T>::maxLegalValue (T value)
00120 {
00121 m_maxLegalValue = value;
00122 }
00123
00124 template <typename T>
00125 inline const T ColumnData<T>::minDataValue () const
00126 {
00127 return m_minDataValue;
00128 }
00129
00130 template <typename T>
00131 inline void ColumnData<T>::minDataValue (T value)
00132 {
00133 m_minDataValue = value;
00134 }
00135
00136 template <typename T>
00137 inline const T ColumnData<T>::maxDataValue () const
00138 {
00139 return m_maxDataValue;
00140 }
00141
00142 template <typename T>
00143 inline void ColumnData<T>::maxDataValue (T value)
00144 {
00145 m_maxDataValue = value;
00146 }
00147
00148 template <typename T>
00149 inline const std::vector<T>& ColumnData<T>::data () const
00150 {
00151 return m_data;
00152 }
00153
00154 template <typename T>
00155 inline void ColumnData<T>::setData (const std::vector<T>& value)
00156 {
00157 m_data = value;
00158 }
00159
00160 template <typename T>
00161 inline T ColumnData<T>::data (int i)
00162 {
00163
00164 return m_data[i - 1];
00165 }
00166
00167 template <typename T>
00168 inline void ColumnData<T>::data (int i, T value)
00169 {
00170
00171 m_data[i - 1] = value;
00172 }
00173
00174
00175
00176 template <typename T>
00177 ColumnData<T>::ColumnData(const ColumnData<T> &right)
00178 :Column(right),
00179 m_minLegalValue(right.m_minLegalValue),
00180 m_maxLegalValue(right.m_maxLegalValue),
00181 m_minDataValue(right.m_minDataValue),
00182 m_maxDataValue(right.m_maxDataValue),
00183 m_data(right.m_data)
00184 {
00185 }
00186
00187 template <typename T>
00188 ColumnData<T>::ColumnData (Table* p)
00189 : Column(p),
00190 m_minLegalValue(),
00191 m_maxLegalValue(),
00192 m_minDataValue(),
00193 m_maxDataValue(),
00194 m_data()
00195 {
00196 }
00197
00198 template <typename T>
00199 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)
00200 : Column(columnIndex,columnName,type,format,unit,p,rpt,w,comment),
00201 m_minLegalValue(),
00202 m_maxLegalValue(),
00203 m_minDataValue(),
00204 m_maxDataValue(),
00205 m_data()
00206 {
00207 }
00208
00209
00210 template <typename T>
00211 ColumnData<T>::~ColumnData()
00212 {
00213 }
00214
00215
00216 template <typename T>
00217 void ColumnData<T>::readColumnData (long firstRow, long nelements, T* nullValue)
00218 {
00219 if ( rows() < nelements )
00220 {
00221 std::cerr << "CCfits: More data requested than contained in table. ";
00222 std::cerr << "Extracting complete column.\n";
00223 nelements = rows();
00224 }
00225
00226 int status(0);
00227 int anynul(0);
00228
00229 FITSUtil::auto_array_ptr<T> array(new T[nelements]);
00230
00231 makeHDUCurrent();
00232
00233 if ( fits_read_col(fitsPointer(),type(), index(), firstRow, 1,
00234 nelements, nullValue, array.get(), &anynul, &status) ) throw FitsError(status);
00235
00236
00237 if (m_data.size() != static_cast<size_t>( rows() ) ) m_data.resize(rows());
00238
00239 std::copy(&array[0],&array[nelements],m_data.begin()+firstRow-1);
00240 if (nelements == rows()) isRead(true);
00241 }
00242
00243 template <typename T>
00244 bool ColumnData<T>::compare (const Column &right) const
00245 {
00246 if ( !Column::compare(right) ) return false;
00247 const ColumnData<T>& that = static_cast<const ColumnData<T>&>(right);
00248 unsigned int n = m_data.size();
00249 if ( that.m_data.size() != n ) return false;
00250 for (unsigned int i = 0; i < n ; i++)
00251 {
00252 if (m_data[i] != that.m_data[i]) return false;
00253 }
00254 return true;
00255 }
00256
00257 template <typename T>
00258 ColumnData<T>* ColumnData<T>::clone () const
00259 {
00260 return new ColumnData<T>(*this);
00261 }
00262
00263 template <typename T>
00264 std::ostream& ColumnData<T>::put (std::ostream& s) const
00265 {
00266 Column::put(s);
00267 if (FITS::verboseMode() && type() != Tstring)
00268 {
00269 s << " Column Legal limits: ( " << m_minLegalValue << "," << m_maxLegalValue << " )\n"
00270 << " Column Data limits: ( " << m_minDataValue << "," << m_maxDataValue << " )\n";
00271 }
00272 if (!m_data.empty())
00273 {
00274 std::ostream_iterator<T> output(s,"\n");
00275
00276
00277 std::copy(m_data.begin(),m_data.end(),output);
00278 }
00279
00280 return s;
00281 }
00282
00283 template <typename T>
00284 void ColumnData<T>::writeData (T* indata, long nRows, long firstRow, T* nullValue)
00285 {
00286
00287
00288
00289
00290
00291 int status(0);
00292 long elementsToWrite(nRows + firstRow -1);
00293
00294 std::vector<T> __tmp(m_data);
00295
00296
00297 if (elementsToWrite != static_cast<long>(m_data.size()))
00298 {
00299
00300 m_data.resize(elementsToWrite,T());
00301 }
00302
00303 std::copy(&indata[0],&indata[nRows],m_data.begin()+firstRow-1);
00304
00305
00306
00307 try
00308 {
00309 if (nullValue)
00310 {
00311 if (fits_write_colnull(fitsPointer(), type(), index(), firstRow, 1, nRows,
00312 indata, nullValue, &status) != 0) throw FitsError(status);
00313 }
00314 else
00315 {
00316 if (fits_write_col(fitsPointer(), type(), index(), firstRow, 1, nRows,
00317 indata, &status) != 0) throw FitsError(status);
00318 }
00319
00320
00321 parent()->updateRows();
00322 }
00323 catch (FitsError)
00324 {
00325
00326 m_data = __tmp;
00327 if (status == NO_NULL) throw NoNullValue(name());
00328 else throw;
00329 }
00330 }
00331
00332 template <typename T>
00333 void ColumnData<T>::writeData (const std::vector<T>& indata, long firstRow, T* nullValue)
00334 {
00335 FITSUtil::CVarray<T> convert;
00336 FITSUtil::auto_array_ptr<T> pcolData (convert(indata));
00337 T* columnData = pcolData.get();
00338 writeData(columnData,indata.size(),firstRow,nullValue);
00339 }
00340
00341 template <typename T>
00342 void ColumnData<T>::insertRows (long first, long number)
00343 {
00344 FITSUtil::FitsNullValue<T> blank;
00345 typename std::vector<T>::iterator in;
00346 if (first !=0)
00347 {
00348 in = m_data.begin()+first;
00349 }
00350 else
00351 {
00352 in = m_data.begin();
00353 }
00354
00355
00356 m_data.insert(in,number,blank());
00357 }
00358
00359 template <typename T>
00360 void ColumnData<T>::deleteRows (long first, long number)
00361 {
00362 m_data.erase(m_data.begin()+first-1,m_data.begin()+first-1+number);
00363 }
00364
00365 template <typename T>
00366 void ColumnData<T>::setDataLimits (T* limits)
00367 {
00368 m_minLegalValue = limits[0];
00369 m_maxLegalValue = limits[1];
00370 m_minDataValue = std::max(limits[2],limits[0]);
00371 m_maxDataValue = std::min(limits[3],limits[1]);
00372 }
00373
00374
00375
00376
00377
00378
00379 #if SPEC_TEMPLATE_IMP_DEFECT || SPEC_TEMPLATE_DECL_DEFECT
00380 template <>
00381 inline void ColumnData<complex<float> >::setDataLimits (complex<float>* limits)
00382 {
00383 m_minLegalValue = limits[0];
00384 m_maxLegalValue = limits[1];
00385 m_minDataValue = limits[2];
00386 m_maxDataValue = limits[3];
00387 }
00388 #else
00389 template <>
00390 void ColumnData<complex<float> >::setDataLimits (complex<float>* limits);
00391 #endif
00392
00393 #if SPEC_TEMPLATE_IMP_DEFECT || SPEC_TEMPLATE_DECL_DEFECT
00394 template <>
00395 inline void ColumnData<complex<double> >::setDataLimits (complex<double>* limits)
00396 {
00397 m_minLegalValue = limits[0];
00398 m_maxLegalValue = limits[1];
00399 m_minDataValue = limits[2];
00400 m_maxDataValue = limits[3];
00401 }
00402 #else
00403 template <>
00404 void ColumnData<complex<double> >::setDataLimits (complex<double>* limits);
00405 #endif
00406
00407
00408 #if SPEC_TEMPLATE_IMP_DEFECT || SPEC_TEMPLATE_DECL_DEFECT
00409 template <>
00410 inline void ColumnData<string>::readColumnData (long firstRow,
00411 long nelements,
00412 string* nullValue)
00413 {
00414 int status = 0;
00415
00416 int anynul = 0;
00417 char** array = new char*[nelements];
00418
00419 int j(0);
00420 for ( ; j < nelements; ++j)
00421 {
00422 array[j] = new char[width() + 1];
00423 }
00424
00425 char* nulval = 0;
00426 if (nullValue)
00427 {
00428 nulval = const_cast<char*>(nullValue->c_str());
00429 }
00430 else
00431 {
00432 nulval = new char;
00433 *nulval = '\0';
00434 }
00435
00436
00437 try
00438 {
00439 makeHDUCurrent();
00440 if (fits_read_col_str(fitsPointer(),index(), firstRow,1,nelements,
00441 nulval,array, &anynul,&status) ) throw FitsError(status);
00442 }
00443 catch (FitsError)
00444 {
00445
00446 for (int jj = 0; jj < nelements; ++jj)
00447 {
00448 delete [] array[jj];
00449 }
00450
00451 delete [] array;
00452 delete nulval;
00453 throw;
00454 }
00455
00456
00457 if (m_data.size() != rows()) setData(std::vector<string>(rows(),string(nulval)));
00458
00459
00460
00461 for ( j = 0; j < nelements; j++)
00462 {
00463 m_data[j - 1 + firstRow] = string(array[j]);
00464 }
00465
00466 for ( j = 0; j < nelements; j++)
00467 {
00468 delete [] array[j];
00469 }
00470
00471 delete [] array;
00472 delete nulval;
00473 if (nelements == rows()) isRead(true);
00474
00475 }
00476 #else
00477 template <>
00478 void ColumnData<string>::readColumnData (long firstRow, long nelements, string* nullValue);
00479 #endif
00480
00481
00482 #if SPEC_TEMPLATE_IMP_DEFECT || SPEC_TEMPLATE_DECL_DEFECT
00483 template <>
00484 inline void ColumnData<complex<float> >::readColumnData (long firstRow,
00485 long nelements,
00486 complex<float>* nullValue)
00487 {
00488
00489 int status(0);
00490 int anynul(0);
00491 FITSUtil::auto_array_ptr<float> pArray(new float[nelements*2]);
00492 float* array = pArray.get();
00493 float nulval(0);
00494 makeHDUCurrent();
00495
00496
00497 if (fits_read_col_cmp(fitsPointer(),index(), firstRow,1,nelements,
00498 nulval,array, &anynul,&status) ) throw FitsError(status);
00499
00500
00501 if (m_data.size() != rows()) m_data.resize(rows());
00502
00503
00504
00505 for (int j = 0; j < nelements; ++j)
00506 {
00507
00508 m_data[j - 1 + firstRow] = std::complex<float>(array[2*j],array[2*j+1]);
00509 }
00510 if (nelements == rows()) isRead(true);
00511
00512 }
00513 #else
00514 template <>
00515 void ColumnData<complex<float> >::readColumnData (long firstRow, long nelements,complex<float>* nullValue );
00516 #endif
00517
00518 #if SPEC_TEMPLATE_IMP_DEFECT || SPEC_TEMPLATE_DECL_DEFECT
00519 template <>
00520 inline void ColumnData<complex<double> >::readColumnData (long firstRow,
00521 long nelements,
00522 complex<double>* nullValue)
00523 {
00524
00525 int status(0);
00526 int anynul(0);
00527 FITSUtil::auto_array_ptr<double> pArray(new double[nelements*2]);
00528 double* array = pArray.get();
00529 double nulval(0);
00530 makeHDUCurrent();
00531
00532
00533 if (fits_read_col_dblcmp(fitsPointer(), index(), firstRow,1,nelements,
00534 nulval,array, &anynul,&status) ) throw FitsError(status);
00535
00536
00537
00538
00539 if (m_data.size() != rows()) setData(std::vector<complex<double> >(rows(),nulval));
00540
00541
00542
00543 for (int j = 0; j < nelements; j++)
00544 {
00545
00546 m_data[j - 1 + firstRow] = std::complex<double>(array[2*j],array[2*j+1]);
00547 }
00548 if (nelements == rows()) isRead(true);
00549
00550 }
00551 #else
00552 template <>
00553 void ColumnData<complex<double> >::readColumnData (long firstRow, long nelements,complex<double>* nullValue);
00554 #endif
00555
00556 #if SPEC_TEMPLATE_DECL_DEFECT
00557 template <>
00558 inline void ColumnData<string>::writeData (const std::vector<string>& indata,
00559 long firstRow, string* nullValue)
00560 {
00561 int status=0;
00562 char** columnData=FITSUtil::CharArray(indata);
00563
00564 if ( fits_write_colnull(fitsPointer(), TSTRING, index(), firstRow, 1, indata.size(),
00565 columnData, 0, &status) != 0 )
00566 throw FitsError(status);
00567 unsigned long elementsToWrite (indata.size() + firstRow - 1);
00568 std::vector<string> __tmp(m_data);
00569 if (m_data.size() < elementsToWrite)
00570 {
00571 m_data.resize(elementsToWrite,"");
00572 std::copy(__tmp.begin(),__tmp.end(),m_data.begin());
00573 }
00574 std::copy(indata.begin(),indata.end(),m_data.begin()+firstRow-1);
00575
00576
00577 for (size_t i = 0; i < indata.size(); ++i)
00578 {
00579 delete [] columnData[i];
00580 }
00581 delete [] columnData;
00582 }
00583 #else
00584 template <>
00585 void ColumnData<string>::writeData (const std::vector<string>& inData, long firstRow, string* nullValue);
00586 #endif
00587
00588 #ifdef SPEC_TEMPLATE_DECL_DEFECT
00589 template <>
00590 inline void ColumnData<complex<float> >::writeData (const std::vector<complex<float> >& inData,
00591 long firstRow,
00592 complex<float>* nullValue)
00593 {
00594 int status(0);
00595 int nRows (inData.size());
00596 FITSUtil::auto_array_ptr<float> pData(new float[nRows*2]);
00597 float* Data = pData.get();
00598 std::vector<complex<float> > __tmp(m_data);
00599 for (int j = 0; j < nRows; ++j)
00600 {
00601 Data[ 2*j] = inData[j].real();
00602 Data[ 2*j + 1] = inData[j].imag();
00603 }
00604
00605 try
00606 {
00607
00608 if (fits_write_col_cmp(fitsPointer(), index(), firstRow, 1,
00609 nRows,Data, &status) != 0) throw FitsError(status);
00610 long elementsToWrite(nRows + firstRow -1);
00611 if (elementsToWrite > static_cast<long>(m_data.size()))
00612 {
00613
00614 m_data.resize(elementsToWrite);
00615 }
00616
00617 std::copy(inData.begin(),inData.end(),m_data.begin()+firstRow-1);
00618
00619
00620 parent()->updateRows();
00621 }
00622 catch (FitsError)
00623 {
00624
00625 m_data.resize(__tmp.size());
00626 m_data = __tmp;
00627 }
00628
00629 }
00630
00631 #else
00632 template <>
00633 void ColumnData<complex<float> >::writeData (const std::vector<complex<float> >& inData, long firstRow,
00634 complex<float>* nullValue);
00635 #endif
00636
00637 #ifdef SPEC_TEMPLATE_DECL_DEFECT
00638 template <>
00639 inline void ColumnData<complex<double> >::writeData (const std::vector<complex<double> >& inData,
00640 long firstRow,
00641 complex<double>* nullValue)
00642 {
00643 int status(0);
00644 int nRows (inData.size());
00645 FITSUtil::auto_array_ptr<double> pData(new double[nRows*2]);
00646 double* Data = pData.get();
00647 std::vector<complex<double> > __tmp(m_data);
00648 for (int j = 0; j < nRows; ++j)
00649 {
00650 pData[ 2*j] = inData[j].real();
00651 pData[ 2*j + 1] = inData[j].imag();
00652 }
00653
00654 try
00655 {
00656
00657 if (fits_write_col_dblcmp(fitsPointer(), index(), firstRow, 1,
00658 nRows,Data, &status) != 0) throw FitsError(status);
00659 long elementsToWrite(nRows + firstRow -1);
00660 if (elementsToWrite > static_cast<long>(m_data.size()))
00661 {
00662
00663 m_data.resize(elementsToWrite);
00664 }
00665
00666 std::copy(inData.begin(),inData.end(),m_data.begin()+firstRow-1);
00667
00668
00669 parent()->updateRows();
00670 }
00671 catch (FitsError)
00672 {
00673
00674 m_data.resize(__tmp.size());
00675 m_data = __tmp;
00676 }
00677
00678 }
00679
00680 #else
00681 template <>
00682 void ColumnData<complex<double> >::writeData (const std::vector<complex<double> >& inData, long firstRow,
00683 complex<double>* nullValue);
00684
00685 #endif
00686 }
00687
00688
00689 #endif