ColumnVectorData.h

00001 //   Read the documentation to learn more about C++ code generator
00002 //   versioning.
00003 //      This is version 1.8 release dated Oct 2007
00004 //      Astrophysics Science Division,
00005 //      NASA/ Goddard Space Flight Center
00006 //      HEASARC
00007 //      http://heasarc.gsfc.nasa.gov
00008 //      e-mail: ccfits@legacy.gsfc.nasa.gov
00009 //
00010 //      Original author: Ben Dorman, L3-Communications EER Systems Inc.
00011 
00012 #ifndef COLUMNVECTORDATA_H
00013 #define COLUMNVECTORDATA_H 1
00014 #ifdef _MSC_VER
00015 #include "MSconfig.h"
00016 #endif
00017 #include "CCfits.h"
00018 
00019 // vector
00020 #include <vector>
00021 // valarray
00022 #include <valarray>
00023 // Column
00024 #include "Column.h"
00025 #ifdef HAVE_CONFIG_H
00026 #include "config.h"
00027 #endif
00028 
00029 #ifdef SSTREAM_DEFECT
00030 #include <strstream>
00031 #else
00032 #include <sstream>
00033 #endif
00034 
00035 #include <memory>
00036 #include <numeric>
00037 namespace CCfits {
00038 
00039         class Table;
00040 
00041 }
00042 
00043 #include "FITS.h"
00044 #include "FITSUtil.h"
00045 using std::complex;
00046 
00047 
00048 namespace CCfits {
00049 
00050 
00051 
00052   template <typename T>
00053   class ColumnVectorData : public Column  //## Inherits: <unnamed>%38BAD1D4D370
00054   {
00055 
00056     public:
00057         ColumnVectorData(const ColumnVectorData< T > &right);
00058         ColumnVectorData (Table* p = 0, T nullVal = FITSUtil::FitsNullValue<T>()());
00059         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 = "", T nullVal = FITSUtil::FitsNullValue<T>()());
00060         ~ColumnVectorData();
00061 
00062         virtual void readData (long firstrow, long nelements, long firstelem = 1);
00063         virtual ColumnVectorData<T>* clone () const;
00064         virtual void setDimen ();
00065         void setDataLimits (T* limits);
00066         const T minLegalValue () const;
00067         void minLegalValue (T value);
00068         const T maxLegalValue () const;
00069         void maxLegalValue (T value);
00070         const T minDataValue () const;
00071         void minDataValue (T value);
00072         const T maxDataValue () const;
00073         void maxDataValue (T value);
00074         const std::vector<std::valarray<T> >& data () const;
00075         void setData (const std::vector<std::valarray<T> >& value);
00076         const std::valarray<T>& data (int i) const;
00077         void data (int i, const std::valarray<T>& value);
00078 
00079       // Additional Public Declarations
00080         friend class Column;
00081     protected:
00082       // Additional Protected Declarations
00083 
00084     private:
00085         ColumnVectorData< T > & operator=(const ColumnVectorData< T > &right);
00086 
00087         virtual bool compare (const Column &right) const;
00088         void resizeDataObject (const std::vector<std::valarray<T> >& indata, size_t firstRow);
00089         //      Reads a specified number of column rows.
00090         //
00091         //      There are no default arguments. The function
00092         //      Column::read(firstrow,firstelem,nelements)
00093         //       is designed for reading the whole column.
00094         virtual void readColumnData (long first, long last, T* nullValue = 0);
00095         virtual std::ostream& put (std::ostream& s) const;
00096         void writeData (const std::valarray<T>& indata, long numRows, long firstRow = 1, T* nullValue = 0);
00097         void writeData (const std::vector<std::valarray<T> >& indata, long firstRow = 1, T* nullValue = 0);
00098         //      Reads a specified number of column rows.
00099         //
00100         //      There are no default arguments. The function
00101         //      Column::read(firstrow,firstelem,nelements)
00102         //       is designed for reading the whole column.
00103         virtual void readRow (size_t row, T* nullValue = 0);
00104         //      Reads a variable row..
00105         virtual void readVariableRow (size_t row, T* nullValue = 0);
00106         void readColumnData (long firstrow, long nelements, long firstelem, T* nullValue = 0);
00107         void writeData (const std::valarray<T>& indata, const std::vector<long>& vectorLengths, long firstRow = 1, T* nullValue = 0);
00108         void writeFixedRow (const std::valarray<T>& data, long row, long firstElem = 1, T* nullValue = 0);
00109         void writeFixedArray (T* data, long nElements, long nRows, long firstRow, T* nullValue = 0);
00110         //      Insert one or more blank rows into a FITS column.
00111         virtual void insertRows (long first, long number = 1);
00112         virtual void deleteRows (long first, long number = 1);
00113         void doWrite (T* array, long row, long rowSize, long firstElem = 1);
00114         const T nullValue () const;
00115         void nullValue (T value);
00116 
00117       // Additional Private Declarations
00118 
00119     private: //## implementation
00120       // Data Members for Class Attributes
00121         T m_nullValue;
00122         T m_minLegalValue;
00123         T m_maxLegalValue;
00124         T m_minDataValue;
00125         T m_maxDataValue;
00126 
00127       // Data Members for Associations
00128         std::vector<std::valarray<T> > m_data;
00129 
00130       // Additional Implementation Declarations
00131 
00132   };
00133 
00134   // Parameterized Class CCfits::ColumnVectorData 
00135 
00136   template <typename T>
00137   inline void ColumnVectorData<T>::readData (long firstrow, long nelements, long firstelem)
00138   {
00139     readColumnData(firstrow,nelements,firstelem,static_cast<T*>(0));
00140   }
00141 
00142   template <typename T>
00143   inline const T ColumnVectorData<T>::nullValue () const
00144   {
00145     return m_nullValue;
00146   }
00147 
00148   template <typename T>
00149   inline void ColumnVectorData<T>::nullValue (T value)
00150   {
00151     m_nullValue = value;
00152   }
00153 
00154   template <typename T>
00155   inline const T ColumnVectorData<T>::minLegalValue () const
00156   {
00157     return m_minLegalValue;
00158   }
00159 
00160   template <typename T>
00161   inline void ColumnVectorData<T>::minLegalValue (T value)
00162   {
00163     m_minLegalValue = value;
00164   }
00165 
00166   template <typename T>
00167   inline const T ColumnVectorData<T>::maxLegalValue () const
00168   {
00169     return m_maxLegalValue;
00170   }
00171 
00172   template <typename T>
00173   inline void ColumnVectorData<T>::maxLegalValue (T value)
00174   {
00175     m_maxLegalValue = value;
00176   }
00177 
00178   template <typename T>
00179   inline const T ColumnVectorData<T>::minDataValue () const
00180   {
00181     return m_minDataValue;
00182   }
00183 
00184   template <typename T>
00185   inline void ColumnVectorData<T>::minDataValue (T value)
00186   {
00187     m_minDataValue = value;
00188   }
00189 
00190   template <typename T>
00191   inline const T ColumnVectorData<T>::maxDataValue () const
00192   {
00193     return m_maxDataValue;
00194   }
00195 
00196   template <typename T>
00197   inline void ColumnVectorData<T>::maxDataValue (T value)
00198   {
00199     m_maxDataValue = value;
00200   }
00201 
00202   template <typename T>
00203   inline const std::vector<std::valarray<T> >& ColumnVectorData<T>::data () const
00204   {
00205     return m_data;
00206   }
00207 
00208   template <typename T>
00209   inline void ColumnVectorData<T>::setData (const std::vector<std::valarray<T> >& value)
00210   {
00211     m_data = value;
00212   }
00213 
00214   template <typename T>
00215   inline const std::valarray<T>& ColumnVectorData<T>::data (int i) const
00216   {
00217     return m_data[i - 1];
00218   }
00219 
00220   template <typename T>
00221   inline void ColumnVectorData<T>::data (int i, const std::valarray<T>& value)
00222   {
00223      if (m_data[i-1].size() != value.size())
00224         m_data[i-1].resize(value.size());
00225      m_data[i - 1] = value;
00226   }
00227 
00228   // Parameterized Class CCfits::ColumnVectorData 
00229 
00230   template <typename T>
00231   ColumnVectorData<T>::ColumnVectorData(const ColumnVectorData<T> &right)
00232       :Column(right),
00233        m_nullValue(right.m_nullValue),
00234        m_minLegalValue(right.m_minLegalValue),
00235        m_maxLegalValue(right.m_maxLegalValue),
00236        m_minDataValue(right.m_minDataValue),
00237        m_maxDataValue(right.m_maxDataValue),
00238        m_data(right.m_data)
00239   {
00240   }
00241 
00242   template <typename T>
00243   ColumnVectorData<T>::ColumnVectorData (Table* p, T nullVal)
00244     : Column(p), m_nullValue(nullVal),
00245        m_minLegalValue(0),
00246        m_maxLegalValue(0),
00247        m_minDataValue(0),
00248        m_maxDataValue(0),
00249        m_data() 
00250   {
00251   }
00252 
00253   template <typename T>
00254   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, T nullVal)
00255         : Column(columnIndex,columnName,type,format,unit,p,rpt,w,comment),
00256           m_nullValue(nullVal),
00257           m_minLegalValue(0),
00258           m_maxLegalValue(0),
00259           m_minDataValue(0),
00260           m_maxDataValue(0), 
00261           m_data()
00262   {
00263   }
00264 
00265 
00266   template <typename T>
00267   ColumnVectorData<T>::~ColumnVectorData()
00268   {
00269   // valarray destructor should do all the work.
00270   }
00271 
00272 
00273   template <typename T>
00274   bool ColumnVectorData<T>::compare (const Column &right) const
00275   {
00276           if ( !Column::compare(right) ) return false;
00277           const ColumnVectorData<T>& that = static_cast<const ColumnVectorData<T>&>(right);
00278           size_t n = m_data.size();
00279           // m_data is of type valarray<T>.
00280           if ( that.m_data.size() != n ) return false;
00281           for (size_t i = 0; i < n ; i++)
00282           {
00283                 size_t nn = m_data[i].size();
00284                 // first check size (also, == on 2 valarrays is only defined if they
00285                 // are equal in size).
00286                 if (that.m_data[i].size() != nn ) return false;
00287 
00288                 std::valarray<bool> test = (m_data[i] == that.m_data[i]);
00289                 for (size_t j = 0; j < nn ; j++ ) if ( !test[j] ) return false;
00290           }
00291           return true;
00292   }
00293 
00294   template <typename T>
00295   ColumnVectorData<T>* ColumnVectorData<T>::clone () const
00296   {
00297   return new ColumnVectorData<T>(*this);
00298   }
00299 
00300   template <typename T>
00301   void ColumnVectorData<T>::resizeDataObject (const std::vector<std::valarray<T> >& indata, size_t firstRow)
00302   {
00303     // the rows() call is the value before updating.
00304     // the updateRows() call at the end sets the call to return the
00305     // value from the fits pointer - which is changed by writeFixedArray
00306     // or writeFixedRow.
00307 
00308     const size_t lastInputRow(indata.size() + firstRow - 1);
00309     const size_t newLastRow = std::max(lastInputRow,static_cast<size_t>(rows()));
00310 
00311     // if the write instruction increases the rows, we need to add
00312     // rows to the data member and preserve its current contents.
00313 
00314     // rows() >= origNRows since it is the value for entire table, 
00315     // not just this column.
00316     const size_t origNRows(m_data.size());
00317     // This will always be an expansion. vector.resize() doesn't
00318     // invalidate any data on an expansion.
00319     if (newLastRow > origNRows) m_data.resize(newLastRow);
00320 
00321     if (varLength())
00322     {
00323        // The incoming data will determine each row size, thus
00324        // no need to preserve any existing values in the row.
00325        // Each value will eventually be overwritten.
00326        for (size_t iRow = firstRow-1; iRow < lastInputRow; ++iRow)
00327        {
00328           std::valarray<T>& current = m_data[iRow];
00329           const size_t newSize = indata[iRow - (firstRow-1)].size();
00330           if (current.size() != newSize)
00331              current.resize(newSize);          
00332        }
00333     }
00334     else
00335     {
00336        // All row sizes in m_data should ALWAYS be either repeat(),
00337        // or 0 if they haven't been initialized.  This is true regardless
00338        // of the incoming data row size.  
00339 
00340        // Perform LAZY initialization of m_data rows.  Only
00341        // expand a row valarray when it is first needed.
00342        for (size_t iRow = firstRow-1; iRow < lastInputRow; ++iRow)
00343        {
00344           if (m_data[iRow].size() != repeat())
00345              m_data[iRow].resize(repeat());
00346        }       
00347     }
00348   }
00349 
00350   template <typename T>
00351   void ColumnVectorData<T>::setDimen ()
00352   {
00353   int status(0);
00354   FITSUtil:: auto_array_ptr<char> dimValue (new char[FLEN_VALUE]);
00355 
00356 #ifdef SSTREAM_DEFECT
00357   std::ostrstream key;
00358 #else
00359   std::ostringstream key;
00360 #endif
00361   key << "TDIM" << index();
00362 
00363 #ifdef SSTREAM_DEFECT
00364   fits_read_key_str(fitsPointer(), key.str(), dimValue.get(),0,&status);
00365 #else
00366   fits_read_key_str(fitsPointer(),const_cast<char*>(key.str().c_str()),dimValue.get(),0,&status);
00367 #endif
00368 
00369   if (status == 0)
00370   {
00371         dimen(String(dimValue.get()));
00372   }
00373   }
00374 
00375   template <typename T>
00376   void ColumnVectorData<T>::readColumnData (long first, long last, T* nullValue)
00377   {
00378   makeHDUCurrent();
00379 
00380 
00381           if ( rows() < last ) 
00382           {
00383                 std::cerr << "CCfits: More data requested than contained in table. ";
00384                 std::cerr << "Extracting complete column.\n";
00385                 last = rows();
00386           }
00387 
00388           long nelements = (last - first + 1)*repeat();
00389 
00390 
00391           readColumnData(first,nelements,1,nullValue);   
00392           if (first <= 1 && last == rows()) isRead(true);
00393   }
00394 
00395   template <typename T>
00396   std::ostream& ColumnVectorData<T>::put (std::ostream& s) const
00397   {
00398   // output header information
00399     Column::put(s);
00400     if ( FITS::verboseMode() )
00401     {
00402           s << " Column Legal limits: ( " << m_minLegalValue << "," << m_maxLegalValue << " )\n" 
00403           << " Column Data  limits: ( " << m_minDataValue << "," << m_maxDataValue << " )\n";
00404     }
00405     if (!m_data.empty())
00406     {
00407           for (size_t j = 0; j < m_data.size(); j++)
00408           {
00409                   size_t n = m_data[j].size();
00410                   if ( n )
00411                   {
00412                           s << "Row " << j + 1 << " Vector Size " << n << '\n';
00413                           for (size_t k = 0; k < n - 1; k++)
00414                           {
00415                                   s << m_data[j][k] << '\t';
00416                           }
00417                           s << m_data[j][n - 1] << '\n';
00418                   }
00419           }
00420     }
00421 
00422     return s;
00423   }
00424 
00425   template <typename T>
00426   void ColumnVectorData<T>::writeData (const std::valarray<T>& indata, long numRows, long firstRow, T* nullValue)
00427   {
00428      // This version of writeData is called by Column write functions which 
00429      // can only write the same number of elements to each row.  
00430      // For fixed width columns, this must be equal to the repeat value
00431      // or an exception is thrown.  For variable width, it only requires
00432      // that indata.size()/numRows is an int.
00433 
00434      // won't do anything if < 0, and will give divide check if == 0.
00435      if (numRows <= 0) throw InvalidNumberOfRows(numRows);
00436 
00437 #ifdef SSTREAM_DEFECT
00438      std::ostrstream msgStr;
00439 #else
00440      std::ostringstream msgStr;
00441 #endif            
00442      if (indata.size() % static_cast<size_t>(numRows))
00443      {
00444         msgStr << "To use this write function, input array size"
00445            <<"\n must be exactly divisible by requested num rows: "
00446            << numRows;
00447         throw InsufficientElements(msgStr.str());
00448      }
00449      const size_t cellsize = indata.size()/static_cast<size_t>(numRows);
00450 
00451      if (!varLength() && cellsize != repeat() )
00452      {      
00453         msgStr << "column: " << name() 
00454                <<  "\n input data size: " << indata.size() 
00455                << " required: " << numRows*repeat();
00456         String msg(msgStr.str());
00457         throw InsufficientElements(msg);     
00458      }
00459 
00460      std::vector<std::valarray<T> > internalFormat(numRows);
00461 
00462      // support writing equal row lengths to variable columns.
00463 
00464      for (long j = 0; j < numRows; ++j)
00465      {
00466         internalFormat[j].resize(cellsize);
00467         internalFormat[j] = indata[std::slice(cellsize*j,cellsize,1)];
00468      }
00469 
00470      // change the size of m_data based on the first row to be written
00471      // and on the input data vector sizes.
00472 
00473      writeData(internalFormat,firstRow,nullValue);    
00474   }
00475 
00476   template <typename T>
00477   void ColumnVectorData<T>::writeData (const std::vector<std::valarray<T> >& indata, long firstRow, T* nullValue)
00478   {
00479      // This is called directly by Column's writeArrays functions, and indirectly
00480      // by both categories of write functions, ie. those which allow differing
00481      // lengths per row and those that don't.
00482     const size_t nInputRows(indata.size());   
00483     using  std::valarray;
00484 
00485     resizeDataObject(indata,firstRow); 
00486     // After the above call, can assume all m_data arrays to be written to 
00487     // have been properly resized whether we're dealing with fixed or
00488     // variable length.       
00489 
00490     if (varLength())
00491     {
00492        // nullValue is ignored for var-length cols, still to keep
00493        // things consistent...
00494        if (nullValue)
00495           m_nullValue = *nullValue;
00496        // firstRow is 1-based, but all these internal row variables 
00497        // will be 0-based.  
00498        const size_t endRow = nInputRows + firstRow-1;
00499        for (size_t iRow = firstRow-1; iRow < endRow; ++iRow)
00500        {
00501           m_data[iRow] = indata[iRow - (firstRow-1)];
00502           // doWrite wants 1-based rows.
00503           doWrite(&m_data[iRow][0], iRow+1, m_data[iRow].size());
00504        }
00505        parent()->updateRows();
00506     }
00507     else
00508     {
00509        // Check for simplest case of all valarrays of size repeat().
00510        // If any are greater, throw an error.
00511        const size_t colRepeat = repeat();
00512        bool allEqualRepeat = true;
00513        for (size_t i=0; i<nInputRows; ++i)
00514        {
00515           const size_t sz = indata[i].size();
00516           if (sz > colRepeat)
00517           {
00518 #ifdef SSTREAM_DEFECT
00519              std::ostrstream oss;
00520 #else
00521              std::ostringstream oss;
00522 #endif 
00523              oss << " vector column length " << colRepeat 
00524                 <<", input valarray length " << sz;
00525              throw InvalidRowParameter(oss.str());               
00526           }
00527           if (sz < colRepeat)
00528              allEqualRepeat = false;
00529        }
00530 
00531        if (allEqualRepeat)
00532        {
00533           // concatenate the valarrays and write.
00534           const size_t nElements (colRepeat*nInputRows);
00535           FITSUtil::CVAarray<T> convert;
00536           FITSUtil::auto_array_ptr<T> pArray(convert(indata));
00537           T* array = pArray.get();
00538 
00539           // if T is complex, then CVAarray returns a 
00540           // C-array of complex objects. But FITS requires an array of complex's
00541           // value_type.
00542 
00543           // This writes to the file and also calls updateRows.
00544           writeFixedArray(array,nElements,nInputRows,firstRow,nullValue);            
00545 
00546           for (size_t j = 0; j < nInputRows ; ++j)
00547           {
00548               const valarray<T>& input   = indata[j];
00549               valarray<T>& current = m_data[j + firstRow - 1];
00550               // current should be resized by resizeDataObject.
00551               current = input;
00552           }
00553        }
00554        else
00555        {
00556           // Some input arrays have fewer than colRepeat elements. 
00557           const size_t endRow = nInputRows + firstRow-1;
00558           for (size_t iRow = firstRow-1; iRow<endRow; ++iRow)
00559           {
00560              // resizeDataObject should already have resized all
00561              // corresponding m_data rows to repeat().
00562              valarray<T>& current = m_data[iRow];
00563              const valarray<T>& input = indata[iRow-(firstRow-1)];
00564              writeFixedRow(input, iRow, 1, nullValue);
00565           }
00566           parent()->updateRows();          
00567        }  
00568 
00569     } // end if !varLength
00570   }
00571 
00572   template <typename T>
00573   void ColumnVectorData<T>::readRow (size_t row, T* nullValue)
00574   {
00575           makeHDUCurrent();
00576 
00577 
00578 
00579           if ( row > static_cast<size_t>(rows()) ) 
00580           {
00581 #ifdef SSTREAM_DEFECT
00582                   std::ostrstream msg;
00583 #else
00584                   std::ostringstream msg;
00585 #endif
00586                 msg << " row requested: " << row << " row range: 1 - " << rows();                
00587 #ifdef SSTREAM_DEFECT
00588                 msg << std::ends;
00589 #endif
00590 
00591                 throw Column::InvalidRowNumber(msg.str()); 
00592           }
00593 
00594           // this is really for documentation purposes. I expect the optimizer will
00595           // remove this redundant definition .
00596           bool variable(type() < 0); 
00597 
00598 
00599           long nelements(repeat());
00600 
00601           if (variable)
00602           {
00603               readVariableRow(row,nullValue);
00604           }
00605           else
00606           {      
00607               readColumnData(row,nelements,1,nullValue);      
00608           }
00609   }
00610 
00611   template <typename T>
00612   void ColumnVectorData<T>::readVariableRow (size_t row, T* nullValue)
00613   {
00614       int status(0);
00615       long offset(0);
00616       long repeat(0);
00617       if (fits_read_descript(fitsPointer(),index(),static_cast<long>(row),
00618                       &repeat,&offset,&status)) throw FitsError(status);
00619       readColumnData(row,repeat,1,nullValue);   
00620   }
00621 
00622   template <typename T>
00623   void ColumnVectorData<T>::readColumnData (long firstrow, long nelements, long firstelem, T* nullValue)
00624   {
00625    int   status=0;
00626 
00627    FITSUtil::auto_array_ptr<T> pArray(new T[nelements]); 
00628    T*     array = pArray.get();
00629    int    anynul(0);
00630 
00631 
00632 
00633    if (fits_read_col(fitsPointer(), abs(type()),index(), firstrow, firstelem,
00634                           nelements, nullValue, array, &anynul, &status) != 0)  
00635        throw FitsError(status);
00636 
00637    size_t countRead = 0;
00638    const size_t ONE = 1;
00639 
00640    if (m_data.size() != static_cast<size_t>(rows())) m_data.resize(rows());
00641    size_t vectorSize(0);
00642    if (!varLength())
00643    {
00644 
00645         vectorSize = std::max(repeat(),ONE); // safety check.
00646 
00647    }
00648    else
00649    {
00650         // assume that the user specified the correct length for 
00651         // variable columns. This should be ok since readVariableColumns
00652         // uses fits_read_descripts to return this information from the
00653         // fits pointer, and this is passed as nelements here.
00654         vectorSize = nelements;       
00655    }
00656    size_t n = nelements; 
00657 
00658    int i = firstrow;
00659    int ii = i - 1;
00660    while ( countRead < n)
00661    {
00662          std::valarray<T>& current = m_data[ii];
00663          if (current.size() != vectorSize) current.resize(vectorSize);
00664          int elementsInFirstRow = vectorSize-firstelem + 1;
00665          bool lastRow = ( (nelements - countRead) < vectorSize);
00666          if (lastRow)
00667          {
00668                int elementsInLastRow = nelements - countRead;
00669                std::valarray<T> ttmp(array + vectorSize*(ii-firstrow) + elementsInFirstRow,
00670                                                      elementsInLastRow);
00671                for (int kk = 0; kk < elementsInLastRow; kk++) current[kk] = ttmp[kk];
00672                countRead += elementsInLastRow;
00673 
00674          }
00675          // what to do with complete rows
00676          else 
00677          {
00678                 if (firstelem == 1 || (firstelem > 1 && i > firstrow) )
00679                 {
00680                         std::valarray<T> ttmp(array + vectorSize*(ii - firstrow) + 
00681                                         elementsInFirstRow,vectorSize);
00682                         current = ttmp;
00683                         ii++;
00684                         i++;
00685                         countRead += vectorSize;   
00686                 }   
00687                 else
00688                 { 
00689                         if (i == firstrow)
00690                         {
00691                                 std::valarray<T> ttmp(array,elementsInFirstRow);
00692                                 for (size_t kk = firstelem ; kk < vectorSize ; kk++)
00693                                       current[kk] = ttmp[kk-firstelem];   
00694                                 countRead += elementsInFirstRow;
00695                                 i++;
00696                                 ii++;
00697                         }
00698                 }
00699          }
00700     }
00701   }
00702 
00703   template <typename T>
00704   void ColumnVectorData<T>::writeData (const std::valarray<T>& indata, const std::vector<long>& vectorLengths, long firstRow, T* nullValue)
00705   {
00706      // Called from Column write functions which allow differing lengths
00707      // for each row.
00708     using namespace std;
00709     const size_t N(vectorLengths.size());
00710     vector<long> sums(N);
00711     // pre-calculate partial sums of vector lengths for use as array offsets.
00712     partial_sum(vectorLengths.begin(),vectorLengths.end(),sums.begin());
00713     // check that sufficient data have been supplied to carry out the entire operation.
00714     if (indata.size() < static_cast<size_t>(sums[N-1]) )
00715     {
00716 #ifdef SSTREAM_DEFECT
00717         ostrstream msgStr;
00718 #else
00719         ostringstream msgStr;
00720 #endif            
00721         msgStr << " input data size: " << indata.size() << " vector length sum: " << sums[N-1];
00722 #ifdef SSTREAM_DEFECT
00723         msgStr << std::ends;
00724 #endif            
00725 
00726         String msg(msgStr.str());
00727         throw InsufficientElements(msg);     
00728     }
00729 
00730     vector<valarray<T> > vvArray(N);
00731     long& last = sums[0];
00732     vvArray[0].resize(last);
00733     for (long jj = 0; jj < last; ++jj) vvArray[0][jj] = indata[jj];
00734 
00735     for (size_t j = 1; j < N; ++j)
00736     {
00737                valarray<T>& __tmp = vvArray[j];
00738                // these  make the code much more readable
00739                long& first = sums[j-1];
00740                long& jlast = sums[j];
00741                __tmp.resize(jlast - first);
00742                for (long k = first; k < jlast; ++k)
00743                { 
00744                         __tmp[k - first] = indata[k];
00745                }
00746     }       
00747 
00748     writeData(vvArray,firstRow,nullValue);
00749   }
00750 
00751   template <typename T>
00752   void ColumnVectorData<T>::writeFixedRow (const std::valarray<T>& data, long row, long firstElem, T* nullValue)
00753   {
00754 
00755     // This is to be called only for FIXED length vector columns.  It will 
00756     // throw if data.size()+firstElem goes beyond the repeat value.
00757     // If data.size() is less than repeat, it leaves the remaining values
00758     // undisturbed both in the file and in m_data storage.
00759 
00760 #ifdef SSTREAM_DEFECT
00761     std::ostrstream msgStr;
00762 #else
00763     std::ostringstream msgStr;
00764 #endif            
00765     if (varLength())
00766     {
00767        msgStr <<"Calling ColumnVectorData::writeFixedRow for a variable length column.\n";
00768        throw FitsFatal(msgStr.str()); 
00769     }
00770 
00771     if (nullValue) m_nullValue = *nullValue;      
00772 
00773     std::valarray<T>& storedRow = m_data[row];    
00774     long inputSize = static_cast<long>(data.size());
00775     long storedSize(storedRow.size());
00776     if (storedSize != static_cast<long>(repeat()))
00777     {
00778        msgStr<<"stored array size vs. column width mismatch in ColumnVectorData::writeFixedRow.\n";
00779        throw FitsFatal(msgStr.str());
00780     }
00781 
00782     if (inputSize + firstElem - 1 > storedSize)
00783     { 
00784           msgStr << " requested write " << firstElem << " to " 
00785                  << firstElem  + inputSize - 1 << " exceeds vector length " << repeat();
00786        throw InvalidRowParameter(msgStr.str());        
00787     }
00788 
00789     // CANNOT give a strong exception safety guarantee because writing
00790     // data changes the file. Any corrective action that could be taken
00791     // [e.g. holding initial contents of the row and writing it back after
00792     // an exception is thrown] could in principle throw the same exception
00793     // we are trying to protect from.
00794 
00795     // routine does however give the weak guarantee (no resource leaks).    
00796 
00797     // It's never a good thing to cast away a const, but doWrite calls the 
00798     // CFITSIO write functions which take a non-const pointer (though
00799     // it shouldn't actually modify the array), and I'd rather not 
00800     // copy the entire valarray just to avoid this problem.
00801     std::valarray<T>& lvData = const_cast<std::valarray<T>&>(data);
00802     T* inPointer = &lvData[0];
00803     doWrite(inPointer, row+1, inputSize, firstElem); 
00804 
00805     // Writing to disk was successful, now update FITS object and return.
00806     const size_t offset = static_cast<size_t>(firstElem) - 1;
00807     for (size_t iElem=0; iElem < static_cast<size_t>(inputSize); ++iElem)
00808     {
00809        // This doesn't require inPointer's non-constness.  It's just
00810        // used here to speed things up a bit.
00811        storedRow[iElem + offset] = inPointer[iElem];
00812     }
00813   }
00814 
00815   template <typename T>
00816   void ColumnVectorData<T>::writeFixedArray (T* data, long nElements, long nRows, long firstRow, T* nullValue)
00817   {
00818     int status(0);
00819 
00820     // check for sanity of inputs, then write to file.
00821     // this function writes only complete rows to a table with
00822     // fixed width rows.
00823 
00824 
00825     if ( nElements < nRows*static_cast<long>(repeat()) )
00826     {
00827 #ifdef SSTREAM_DEFECT
00828         std::ostrstream msgStr;
00829 #else
00830         std::ostringstream msgStr;
00831 #endif
00832         msgStr << " input array size: " << nElements << " required " << nRows*repeat();
00833         String msg(msgStr.str());
00834 
00835             throw Column::InsufficientElements(msg);
00836     } 
00837 
00838     if (nullValue) m_nullValue = *nullValue;      
00839 
00840     if (fits_write_colnull(fitsPointer(),abs(type()),index(),firstRow,
00841                         1,nElements,data,nullValue,&status)) throw FitsError(status);
00842 
00843     parent()->updateRows();
00844   }
00845 
00846   template <typename T>
00847   void ColumnVectorData<T>::insertRows (long first, long number)
00848   {
00849     typename std::vector<std::valarray<T> >::iterator in;
00850     if (first !=0) 
00851     {
00852             in = m_data.begin()+first;
00853     }
00854     else
00855     {
00856             in = m_data.begin();
00857     }           
00858 
00859     // non-throwing operations.
00860     m_data.insert(in,number,std::valarray<T>(T(),0));
00861   }
00862 
00863   template <typename T>
00864   void ColumnVectorData<T>::deleteRows (long first, long number)
00865   {
00866     // the following is an ugly workaround for a bug in g++ v3.0 that
00867     // does not erase vector elements cleanly in this case.
00868 
00869     size_t N(m_data.size());
00870     int newSize(N - number);      
00871     std::vector<std::valarray<T> > __tmp(newSize);
00872 
00873     int lastDeleted( number + first - 1 );
00874     int firstDeleted(first);
00875     int count(0);
00876         {
00877                 for ( size_t j = 1; j <= N; ++j)
00878                 {
00879                         if (  (j - firstDeleted)*(lastDeleted - j) >= 0 )       
00880                         {                ++count; 
00881                         } 
00882                         else
00883                         {
00884                                 __tmp[j - 1 - count].resize(m_data[j - 1].size());
00885                                 __tmp[j - 1 - count] = m_data[j - 1];
00886                         }
00887                 }                           
00888     }
00889 
00890     m_data.clear();
00891     m_data.resize(newSize);
00892         {
00893                 for (int j = 0; j < newSize; ++j)
00894                 {
00895                         m_data[j].resize(__tmp[j].size());
00896                         m_data[j] = __tmp[j];
00897                 }
00898         }
00899   }
00900 
00901   template <typename T>
00902   void ColumnVectorData<T>::setDataLimits (T* limits)
00903   {
00904     m_minLegalValue = limits[0];
00905     m_maxLegalValue = limits[1];
00906     m_minDataValue = std::max(limits[2],limits[0]);
00907     m_maxDataValue = std::min(limits[3],limits[1]);
00908   }
00909 
00910   template <typename T>
00911   void ColumnVectorData<T>::doWrite (T* array, long row, long rowSize, long firstElem)
00912   {
00913     int status(0);
00914     // internal functioning of write_colnull forbids its use for writing
00915     // variable width columns. If a nullvalue argument was supplied it will
00916     // be ignored.
00917     if ( !varLength())
00918     {
00919         if (fits_write_colnull(fitsPointer(),type(),index(),row, firstElem, rowSize,
00920                     array,&m_nullValue,&status)) throw FitsError(status);
00921     }
00922     else
00923     {
00924         if (fits_write_col(fitsPointer(),abs(type()),index(),row,firstElem,rowSize,
00925                     array,&status)) throw FitsError(status);
00926 
00927     }
00928   }
00929 
00930   // Additional Declarations
00931 
00932   // all functions that operate on complex data that call cfitsio 
00933   // need to be specialized. The signature containing complex<T>* objects
00934   // is unfortunate, perhaps, for this purpose, but the user will  access
00935   // rw operations through standard library containers.
00936 
00937 
00938 
00939 
00940 
00941 #if SPEC_TEMPLATE_IMP_DEFECT || SPEC_TEMPLATE_DECL_DEFECT
00942 template <>
00943 inline void ColumnVectorData<complex<float> >::setDataLimits (complex<float>* limits)
00944         {
00945                 m_minLegalValue = limits[0];
00946                 m_maxLegalValue = limits[1];
00947                 m_minDataValue =  limits[2];
00948                 m_maxDataValue =  limits[3];
00949         }
00950 #else
00951 template <>
00952   void 
00953   ColumnVectorData<complex<float> >::setDataLimits (complex<float>* limits);
00954 #endif
00955 
00956 #if SPEC_TEMPLATE_IMP_DEFECT || SPEC_TEMPLATE_DECL_DEFECT
00957 template <>
00958 inline void ColumnVectorData<complex<double> >::setDataLimits (complex<double>* limits)
00959         {
00960                 m_minLegalValue = limits[0];
00961                 m_maxLegalValue = limits[1];
00962                 m_minDataValue =  limits[2];
00963                 m_maxDataValue =  limits[3];
00964         }
00965 #else
00966  template <>
00967    void 
00968    ColumnVectorData<complex<double> >::setDataLimits (complex<double>* limits);
00969 #endif
00970 
00971 
00972 #if SPEC_TEMPLATE_IMP_DEFECT || SPEC_TEMPLATE_DECL_DEFECT
00973         template <>
00974         inline void ColumnVectorData<std::complex<float> >::readColumnData(long firstRow, 
00975                                 long nelements, long firstElem, std::complex<float>* null )
00976         {
00977             int   status=0;
00978             float nulval (0);
00979             FITSUtil::auto_array_ptr<float> pArray(new float[2*nelements]); 
00980             float*     array = pArray.get();
00981             int    anynul(0);
00982 
00983             if (fits_read_col_cmp(fitsPointer(),index(),firstRow, firstElem,
00984                             nelements,nulval,array,&anynul,&status) ) throw FitsError(status);
00985 
00986             if (m_data.size() != static_cast<size_t>(rows())) m_data.resize(rows());
00987 
00988             std::valarray<std::complex<float> > readData(nelements);
00989             for (long j = 0; j < nelements; ++j)
00990             {
00991                     readData[j] = std::complex<float>(array[2*j],array[2*j+1]);
00992             }
00993             size_t countRead = 0;
00994             const size_t ONE = 1;
00995 
00996             if (m_data.size() != static_cast<size_t>(rows())) m_data.resize(rows());
00997             size_t vectorSize(0);
00998             if (!varLength())
00999             {
01000                  vectorSize = std::max(repeat(),ONE); // safety check.
01001             }
01002             else
01003             {
01004                  // assume that the user specified the correct length for 
01005                  // variable columns. This should be ok since readVariableColumns
01006                  // uses fits_read_descripts to return this information from the
01007                  // fits pointer, and this is passed as nelements here.
01008                  vectorSize = nelements;       
01009             }
01010             size_t n = nelements; 
01011 
01012             int i = firstRow;
01013             int ii = i - 1;
01014             while ( countRead < n)
01015             {
01016                     std::valarray<complex<float> >& current = m_data[ii];
01017                     if (current.size() != vectorSize) current.resize(vectorSize,0.);
01018                     int elementsInFirstRow = vectorSize-firstElem + 1;
01019                     bool lastRow = ( (nelements - countRead) < vectorSize);
01020                     if (lastRow)
01021                     {
01022                             int elementsInLastRow = nelements - countRead;
01023                             std::copy(&readData[countRead],&readData[nelements],&current[0]);
01024                             countRead += elementsInLastRow;
01025                     }             
01026                     // what to do with complete rows. if firstElem == 1 the 
01027                     else 
01028                     {
01029                             if (firstElem == 1 || (firstElem > 1 && i > firstRow) )
01030                             {
01031                                     current = readData[std::slice(vectorSize*(ii-firstRow)+
01032                                                                elementsInFirstRow,vectorSize,1)];
01033                                     ++ii;
01034                                     ++i;
01035                                     countRead += vectorSize;   
01036                             }   
01037                             else
01038                             { 
01039                                     if (i == firstRow)
01040                                     {
01041                                             std::copy(&readData[0],&readData[elementsInFirstRow],
01042                                                                             &current[firstElem]);
01043                                             countRead += elementsInFirstRow;
01044                                             ++i;
01045                                             ++ii;
01046                                     }
01047                             }
01048                     }
01049             }
01050     }
01051 #else
01052 template <>
01053 void ColumnVectorData<complex<float> >::readColumnData(long firstRow, 
01054                         long nelements, 
01055                         long firstElem, complex<float>* null);
01056 #endif
01057 
01058 #if SPEC_TEMPLATE_IMP_DEFECT || SPEC_TEMPLATE_DECL_DEFECT
01059     template <>
01060     inline void ColumnVectorData<complex<double> >::readColumnData (long firstRow, 
01061               long nelements,long firstElem, 
01062               complex<double>* nullValue)
01063     {
01064 
01065         // duplicated for each complex type to work around imagined or
01066         // actual compiler deficiencies.
01067             int   status=0;
01068             double nulval (0);
01069             FITSUtil::auto_array_ptr<double> pArray(new double[2*nelements]); 
01070             double*     array = pArray.get();
01071             int    anynul(0);
01072 
01073             if (fits_read_col_dblcmp(fitsPointer(),index(),firstRow, firstElem,
01074                             nelements,nulval,array,&anynul,&status) ) throw FitsError(status);
01075 
01076             if (m_data.size() != static_cast<size_t>(rows())) m_data.resize(rows());
01077 
01078             std::valarray<std::complex<double> > readData(nelements);
01079             for (long j = 0; j < nelements; ++j)
01080             {
01081                     readData[j] = std::complex<double>(array[2*j],array[2*j+1]);
01082             }
01083             size_t countRead = 0;
01084             const size_t ONE = 1;
01085 
01086             if (m_data.size() != static_cast<size_t>(rows())) m_data.resize(rows());
01087             size_t vectorSize(0);
01088             if (!varLength())
01089             {
01090                  vectorSize = std::max(repeat(),ONE); // safety check.
01091             }
01092             else
01093             {
01094                  // assume that the user specified the correct length for 
01095                  // variable columns. This should be ok since readVariableColumns
01096                  // uses fits_read_descripts to return this information from the
01097                  // fits pointer, and this is passed as nelements here.
01098                  vectorSize = nelements;       
01099             }
01100             size_t n = nelements; 
01101 
01102             int i = firstRow;
01103             int ii = i - 1;
01104             while ( countRead < n)
01105             {
01106                     std::valarray<std::complex<double> >& current = m_data[ii];
01107                     if (current.size() != vectorSize) current.resize(vectorSize,0.);
01108                     int elementsInFirstRow = vectorSize-firstElem + 1;
01109                     bool lastRow = ( (nelements - countRead) < vectorSize);
01110                     if (lastRow)
01111                     {
01112                             int elementsInLastRow = nelements - countRead;
01113                             std::copy(&readData[countRead],&readData[nelements],&current[0]);
01114                             countRead += elementsInLastRow;
01115                     }             
01116                     // what to do with complete rows. if firstElem == 1 the 
01117                     else 
01118                     {
01119                             if (firstElem == 1 || (firstElem > 1 && i > firstRow) )
01120                             {
01121                                     current = readData[std::slice(vectorSize*(ii-firstRow)+
01122                                                                elementsInFirstRow,vectorSize,1)];
01123                                     ++ii;
01124                                     ++i;
01125                                     countRead += vectorSize;   
01126                             }   
01127                             else
01128                             { 
01129                                     if (i == firstRow)
01130                                     {
01131                                             std::copy(&readData[0],&readData[elementsInFirstRow],
01132                                                                             &current[firstElem]);
01133                                             countRead += elementsInFirstRow;
01134                                             ++i;
01135                                             ++ii;
01136                                     }
01137                             }
01138                     }
01139             }
01140     }
01141 #else
01142 template <>
01143 void ColumnVectorData<complex<double> >::readColumnData (long firstRow, 
01144                         long nelements,
01145                         long firstElem, complex<double>* null);
01146 #endif
01147 
01148 #if SPEC_TEMPLATE_IMP_DEFECT || SPEC_TEMPLATE_DECL_DEFECT
01149         template <>
01150         inline void ColumnVectorData<complex<float> >::writeFixedArray 
01151                         (complex<float>* data, long nElements, long nRows, long firstRow, 
01152                          complex<float>* nullValue)
01153         {
01154 
01155                 int status(0);
01156 
01157     // check for sanity of inputs, then write to file.
01158     // this function writes only complete rows to a table with
01159     // fixed width rows.
01160 
01161 
01162                 if ( nElements < nRows*static_cast<long>(repeat()) )
01163                 {
01164 #ifdef SSTREAM_DEFECT
01165                         std::ostrstream msgStr;
01166 #else
01167                         std::ostringstream msgStr;
01168 #endif
01169                         msgStr << " input array size: " << nElements 
01170                                         << " required " << nRows*repeat();
01171 #ifdef SSTREAM_DEFECT
01172                         msgStr << std::ends;
01173 #endif
01174 
01175 
01176                         String msg(msgStr.str());
01177 
01178                         throw Column::InsufficientElements(msg);
01179                 } 
01180 
01181                 FITSUtil::auto_array_ptr<float> realData(new float[2*nElements]);
01182 
01183                 for (int j = 0; j < nElements; ++j)
01184                 {
01185                         realData[2*j] = data[j].real();
01186                         realData[2*j+1] = data[j].imag();       
01187                 }
01188 
01189 
01190 
01191                 if (fits_write_col_cmp(fitsPointer(),index(),firstRow,
01192                         1,nElements,realData.get(),&status)) throw FitsError(status);
01193 
01194                 parent()->updateRows();
01195         }
01196 #else
01197 template <>
01198 void ColumnVectorData<complex<float> >::writeFixedArray 
01199      (complex<float>* data, long nElements, long nRows, long firstRow, std::complex<float>* null);
01200 #endif
01201 
01202 #if SPEC_TEMPLATE_IMP_DEFECT || SPEC_TEMPLATE_DECL_DEFECT
01203         template <>
01204         inline void ColumnVectorData<complex<double> >::writeFixedArray 
01205                         (complex<double>* data, long nElements, long nRows, long firstRow, 
01206                          complex<double>* nullValue)
01207         {
01208                 int status(0);
01209 
01210     // check for sanity of inputs, then write to file.
01211     // this function writes only complete rows to a table with
01212     // fixed width rows.
01213 
01214 
01215                 if ( nElements < nRows*static_cast<long>(repeat()) )
01216                 {
01217 #ifdef SSTREAM_DEFECT
01218                         std::ostrstream msgStr;
01219 #else
01220                         std::ostringstream msgStr;
01221 #endif
01222                         msgStr << " input array size: " << nElements 
01223                                         << " required " << nRows*repeat();
01224 #ifdef SSTREAM_DEFECT
01225                         msgStr << std::ends;
01226 #endif
01227 
01228                         String msg(msgStr.str());
01229 
01230                         throw Column::InsufficientElements(msg);
01231                 } 
01232 
01233                 FITSUtil::auto_array_ptr<double> realData(new double[2*nElements]);
01234 
01235                 for (int j = 0; j < nElements; ++j)
01236                 {
01237                         realData[2*j] = data[j].real();
01238                         realData[2*j+1] = data[j].imag();       
01239                 }
01240 
01241 
01242 
01243                 if (fits_write_col_dblcmp(fitsPointer(),index(),firstRow,
01244                         1,nElements,realData.get(),&status)) throw FitsError(status);
01245 
01246                 parent()->updateRows();
01247 
01248         }
01249 #else
01250 template <>
01251 void ColumnVectorData<complex<double> >::writeFixedArray 
01252                 (complex<double>* data, long nElements, long nRows, long firstRow, 
01253                  std::complex<double>* null);
01254 #endif
01255 
01256 #ifdef SPEC_TEMPLATE_DECL_DEFECT
01257   template <>
01258   inline void  
01259   ColumnVectorData<std::complex<float> >::doWrite 
01260   (std::complex<float>* data, long row, long rowSize, long firstElem )
01261   {
01262     int status(0);
01263     FITSUtil::auto_array_ptr<float> carray( new float[2*rowSize]); 
01264     for ( long j = 0 ; j < rowSize; ++ j)
01265       {
01266         carray[2*j] = data[j].real();
01267         carray[2*j + 1] = data[j].imag();
01268       }
01269     if (fits_write_col_cmp(fitsPointer(),index(),row,firstElem,rowSize,
01270                            carray.get(),&status)) throw FitsError(status);
01271   }
01272 
01273 
01274   template <>
01275   inline void  
01276   ColumnVectorData<std::complex<double> >::doWrite
01277   (std::complex<double>* data, long row, long rowSize, long firstElem )
01278   {
01279     int status(0);
01280     FITSUtil::auto_array_ptr<double> carray( new double[2*rowSize]); 
01281     for ( long j = 0 ; j < rowSize; ++ j)
01282       {
01283         carray[2*j] = data[j].real();
01284         carray[2*j + 1] = data[j].imag();
01285       }
01286     if (fits_write_col_dblcmp(fitsPointer(),index(),row,firstElem,rowSize,
01287                               carray.get(),&status)) throw FitsError(status);
01288 
01289   }
01290 
01291 #else
01292 template<>
01293 void 
01294 ColumnVectorData<complex<float> >::doWrite 
01295                 ( complex<float>* data, long row, long rowSize, long firstElem );
01296 
01297 template<>
01298 void 
01299 ColumnVectorData<complex<double> >::doWrite 
01300                 ( complex<double>* data, long row, long rowSize, long firstElem );
01301 #endif
01302 } // namespace CCfits
01303 
01304 
01305 #endif

Generated on Fri Oct 12 13:39:38 2007 for CCfits by  doxygen 1.4.7