Index: /trunk/Mars/mcore/fits.h =================================================================== --- /trunk/Mars/mcore/fits.h (revision 16283) +++ /trunk/Mars/mcore/fits.h (revision 16284) @@ -69,4 +69,12 @@ { public: + //I know I know, you're going to yiell that this does not belong here. + //It will belong in the global scope eventually, and it makes the coding of zfits much simpler this way. + typedef enum + { + UNCOMPRESSED, + SMOOTHMAN + } FitsCompression; + struct Entry { @@ -74,4 +82,5 @@ string value; string comment; + string fitsString; template @@ -90,4 +99,6 @@ { off_t offset; + + bool isCompressed; string name; @@ -103,10 +114,13 @@ char type; string unit; + FitsCompression comp; }; typedef map Keys; typedef map Columns; + typedef vector SortedColumns; Columns cols; + SortedColumns sortedCols; Keys keys; @@ -200,5 +214,14 @@ // Set value, comment and type - com = ppp==string::npos ? "" : Trim(val.substr(ppp+1)); + if (ppp==string::npos) + com = ""; + else + {// comments could be just spaces. take care of this. + if (val.size() == ppp+1) + com = ""; + else + com = Trim(val.substr(ppp+1)); + } +// com = ppp==string::npos ? "" : Trim(val.substr(ppp+1)); val = Trim(val.substr(1, p-2)); type = 'T'; @@ -208,5 +231,8 @@ const size_t p = val.find_first_of('/'); - com = Trim(val.substr(p+2)); + if (val.size() == p+1) + com = ""; + else + com = Trim(val.substr(p+2)); val = Trim(val.substr(0, p)); @@ -217,5 +243,5 @@ } - const Entry e = { type, val, com }; + const Entry e = { type, val, com, vec[i] }; rc[key] = e; } @@ -226,10 +252,15 @@ Table() : offset(0) { } Table(const vector &vec, off_t off) : - offset(off), keys(ParseBlock(vec)) - { + offset(off), isCompressed(false), keys(ParseBlock(vec)) + { + if (HasKey("ZTABLE")) + if (Check("ZTABLE", 'B', "T")) + { + isCompressed = true; + } + if (!Check("XTENSION", 'T', "BINTABLE") || !Check("NAXIS", 'I', "2") || !Check("BITPIX", 'I', "8") || - !Check("PCOUNT", 'I', "0") || !Check("GCOUNT", 'I', "1") || !Check("EXTNAME", 'T') || @@ -239,10 +270,22 @@ return; - bytes_per_row = Get("NAXIS1"); - num_rows = Get("NAXIS2"); + if (isCompressed) + { + if (!Check("ZPCOUNT", 'I', "0")) + return; + } + else + { + if (!Check("PCOUNT", 'I', "0")) + return; + } + + bytes_per_row = isCompressed ? Get("ZNAXIS1") : Get("NAXIS1"); + num_rows = isCompressed ? Get("ZNAXIS2") : Get("NAXIS2"); num_cols = Get("TFIELDS"); - datasum = Get("DATASUM", -1); + datasum = isCompressed ? Get("ZDATASUM", -1) : Get("DATASUM", -1); size_t bytes = 0; + string tFormName = isCompressed ? "ZFORM" : "TFORM"; for (size_t i=1; i<=num_cols; i++) { @@ -251,10 +294,15 @@ if (!Check("TTYPE"+num.str(), 'T') || - !Check("TFORM"+num.str(), 'T')) + !Check(tFormName+num.str(), 'T')) return; const string id = Get("TTYPE"+num.str()); - const string fmt = Get("TFORM"+num.str()); + const string fmt = Get(tFormName+num.str()); const string unit = Get("TUNIT"+num.str(), ""); + const string comp = Get("ZCTYP"+num.str(), ""); + + FitsCompression compress = UNCOMPRESSED; + if (comp == "SMOO") + compress = SMOOTHMAN; istringstream sin(fmt); @@ -297,7 +345,8 @@ } - const Table::Column col = { bytes, n, size, type, unit }; + const Table::Column col = { bytes, n, size, type, unit, compress}; cols[id] = col; + sortedCols.push_back(col); bytes += n*size; } @@ -420,5 +469,5 @@ }; -private: +protected: ofstream fCopy; @@ -500,5 +549,46 @@ } - void Constructor(const string &fname, string fout, bool force) + //There may be a gap between the main table and the start of the heap: this computes the offset + size_t ComputeGapFromDataToHeap(const Table& table) + { + //get row width + const size_t rowWidth = table.Get("NAXIS1"); + const size_t numRows = table.Get("NAXIS2"); + + //get offset of special data area from start of main table + size_t heapShift = HasKey("THEAP") ? table.Get("THEAP") : 0; + if (heapShift > rowWidth*numRows) + heapShift -= rowWidth*numRows; + else + heapShift = 0; + + return heapShift; + } + + //skip the current table + bool SkipCurrentTable(const vector& vector) + { + Table currentTable(vector, tellg()); + //get row width + const size_t rowWidth = currentTable.Get("NAXIS1"); + const size_t numRows = currentTable.Get("NAXIS2"); + + //get offset of special data area from start of main table + const off_t heapShift = ComputeGapFromDataToHeap(currentTable); + + //and special data area size + const off_t heapSize = HasKey("PCOUNT") ? currentTable.Get("PCOUNT") : 0; + + //skip the row of the table, padded to 2880 bytes + off_t paddingSize = (rowWidth*numRows + heapSize + heapShift)%2880; + paddingSize += (paddingSize==0) ? 0 : 2880 - paddingSize; + + const off_t whereDoIGo = tellg() + (off_t)(rowWidth*numRows) + heapSize + heapShift + paddingSize; + seekg(whereDoIGo); + + return good(); + } + + void Constructor(const string &fname, string fout, const string& tableName, bool force, bool mightFail) { char simple[10]; @@ -529,4 +619,7 @@ if (!good()) { + //we allow that the table is not found (to be checked later on with bool casting) + if (mightFail) { return;} + clear(rdstate()|ios::badbit); #ifdef __EXCEPTIONS @@ -565,6 +658,16 @@ if (block[0].substr(0, 9)=="XTENSION=") { - // FIXME: Check for table name - + //Check for table name. Skip until eof or requested table are found. + if (tableName != "") + { + Table currentTable(block, tellg()); + string currentName = currentTable.Get("EXTNAME"); + if (currentName != tableName) + { + SkipCurrentTable(block); + fChkHeader.reset(); + continue; + } + } fTable = Table(block, tellg()); fRow = (size_t)-1; @@ -630,12 +733,12 @@ public: - fits(const string &fname, bool force=false) : izstream(fname.c_str()) - { - Constructor(fname, "", force); - } - - fits(const string &fname, const string &fout, bool force=false) : izstream(fname.c_str()) - { - Constructor(fname, fout, force); + fits(const string &fname, const string& tableName="", bool force=false, bool mightFail=false) : izstream(fname.c_str()) + { + Constructor(fname, "", tableName, force, mightFail); + } + + fits(const string &fname, const string &fout, const string& tableName="", bool force=false, bool mightFail=false) : izstream(fname.c_str()) + { + Constructor(fname, fout, tableName, force, mightFail); } @@ -647,9 +750,14 @@ } - uint8_t ReadRow(size_t row) + virtual void StageRow(size_t row, char* dest) { // if (row!=fRow+1) // Fast seeking is ensured by izstream seekg(fTable.offset+row*fTable.bytes_per_row); - + read(dest, fTable.bytes_per_row); + //fin.clear(fin.rdstate()&~ios::eofbit); + } + + uint8_t ReadRow(size_t row) + { // For the checksum we need everything to be correctly aligned const uint8_t offset = (row*fTable.bytes_per_row)%4; @@ -671,6 +779,5 @@ *--ie = 0; - read(fBufferRow.data()+offset, fTable.bytes_per_row); - //fin.clear(fin.rdstate()&~ios::eofbit); + StageRow(row, fBufferRow.data()+offset); if (row==fRow+1) @@ -699,8 +806,21 @@ } + virtual void MoveColumnDataToUserSpace(char* dest, const char*src, const Table::Column& c) + { + // Let the compiler do some optimization by + // knowing that we only have 1, 2, 4 and 8 + switch (c.size) + { + case 1: memcpy (dest, src, c.num*c.size); break; + case 2: revcpy<2>(dest, src, c.num); break; + case 4: revcpy<4>(dest, src, c.num); break; + case 8: revcpy<8>(dest, src, c.num); break; + } + } + #if !defined(__MARS__) && !defined(__CINT__) - bool GetRow(size_t row, bool check=true) -#else - bool GetRowNum(size_t row, bool check=true) + virtual bool GetRow(size_t row, bool check=true) +#else + virtual bool GetRowNum(size_t row, bool check=true) #endif { @@ -721,13 +841,5 @@ char *dest = reinterpret_cast(it->first); - // Let the compiler do some optimization by - // knowing that we only have 1, 2, 4 and 8 - switch (c.size) - { - case 1: memcpy (dest, src, c.num*c.size); break; - case 2: revcpy<2>(dest, src, c.num); break; - case 4: revcpy<4>(dest, src, c.num); break; - case 8: revcpy<8>(dest, src, c.num); break; - } + MoveColumnDataToUserSpace(dest, src, c); } @@ -744,5 +856,5 @@ } - bool SkipNextRow() + virtual bool SkipNextRow() { seekg(fTable.offset+(++fRow)*fTable.bytes_per_row); @@ -891,4 +1003,5 @@ bool HasColumn(const string& col) const { return fTable.HasColumn(col);} const Table::Columns &GetColumns() const { return fTable.GetColumns();} + const Table::SortedColumns& GetSortedColumns() const { return fTable.sortedCols;} const Table::Keys &GetKeys() const { return fTable.GetKeys();} @@ -914,4 +1027,5 @@ bool IsFileOk() const { return (fChkHeader+fChkData).valid(); } + bool IsCompressedFITS() const { return fTable.isCompressed;} };