Changeset 16814 for trunk/Mars/mcore


Ignore:
Timestamp:
06/12/13 15:56:12 (11 years ago)
Author:
lyard
Message:
added block headers and changed compression mechanism
File:
1 edited

Legend:

Unmodified
Added
Removed
  • trunk/Mars/mcore/zfits.h

    r16811 r16814  
    1313#include "fits.h"
    1414#include "huffman.h"
     15
     16
     17#define FACT_RAW       0x0
     18#define FACT_SMOOTHING 0x1
     19#define FACT_HUFFMAN16 0x2
     20
     21#define FACT_COL_MAJOR 'C'
     22#define FACT_ROW_MAJOR 'R'
    1523
    1624
     
    7987
    8088private:
     89
     90    //Structure helper for reading tiles headers
    8191    typedef struct TileHeader
    8292    {
     
    8494      uint32_t numRows;
    8595      uint64_t size;
    86       TileHeader(uint32_t nRows=0,
    87                  uint64_t s=0) : id({'T', 'I', 'L', 'E'}),
     96
     97      TileHeader() {}
     98
     99      TileHeader(uint32_t nRows,
     100                 uint64_t s) : id({'T', 'I', 'L', 'E'}),
    88101                                 numRows(nRows),
    89102                                 size(s)
    90103      { };
    91       friend ostream& operator << (ostream& out, const TileHeader& h)
    92       {
    93           out << h.id[0] << h.id[1] << h.id[2] << h.id[3] << " num Rows: " << h.numRows << ", tile size: " << h.size;
    94           return out;
    95       }
    96104    } __attribute__((__packed__)) TileHeader;
     105
     106    //Structure helper for reading blocks headers and compresion schemes
     107    typedef struct BlockHeader
     108    {
     109        uint64_t      size;
     110        char          ordering;
     111        unsigned char numProcs;
     112        uint16_t      processings[];
     113
     114        BlockHeader(uint64_t      s=0,
     115                    char          o=FACT_ROW_MAJOR,
     116                    unsigned char n=1) : size(s),
     117                                         ordering(o),
     118                                         numProcs(n)
     119        {}
     120    } __attribute__((__packed__)) BlockHeader;
    97121
    98122    // Do what it takes to initialize the compressed structured
     
    103127            return;
    104128
     129        if (fTable.isCompressed)
     130        for (auto it=fTable.sortedCols.begin(); it!= fTable.sortedCols.end(); it++)
     131            if (it->comp != FACT)
     132            {
     133#ifdef __EXCEPTIONS
     134                    throw runtime_error("ERROR: Only the FACT compression scheme is handled by this reader.");
     135#else
     136                    gLog << ___err___ << "ERROR: Only the FACT compression scheme is handled by this reader." << endl;
     137                    return;
     138#endif
     139            }
     140
     141        fColumnOrdering.resize(fTable.sortedCols.size());
     142        for (auto it=fColumnOrdering.begin(); it != fColumnOrdering.end(); it++)
     143            (*it) = FACT_ROW_MAJOR;
    105144        //Get compressed specific keywords
    106145        fNumTiles       = fTable.isCompressed ? GetInt("NAXIS2") : 0;
     
    129168    vector<char> fTransposedBuffer; ///<intermediate buffer to transpose the rows
    130169    vector<char> fCompressedBuffer; ///<compressed rows
     170    vector<char> fColumnOrdering; ///< ordering of the column's rows
    131171
    132172    size_t fNumTiles;       ///< Total number of tiles
     
    137177
    138178    vector<vector<pair<int64_t, int64_t> > > fCatalog;///< Catalog, i.e. the main table that points to the compressed data.
    139     vector<size_t> fTileSize; ///< size in bytes of each compressed tile
    140     vector<vector<size_t> > fTileOffsets; ///< offset from start of tile of a given compressed column
     179    vector<size_t>                           fTileSize; ///< size in bytes of each compressed tile
     180    vector<vector<size_t> >                  fTileOffsets; ///< offset from start of tile of a given compressed column
    141181
    142182    void AllocateBuffers()
     
    148188
    149189        fTransposedBuffer.resize(fTable.bytes_per_row*fNumRowsPerTile);
    150         fCompressedBuffer.resize(fTable.bytes_per_row*fNumRowsPerTile + fTable.num_cols); //use a bit more memory for compression flags
     190        fCompressedBuffer.resize(fTable.bytes_per_row*fNumRowsPerTile + fTable.num_cols*(sizeof(BlockHeader)+256)); //use a bit more memory for block headers
    151191    }
    152192
     
    266306            const char *src = fTransposedBuffer.data();
    267307
    268             for (auto it=fTable.sortedCols.begin(); it!=fTable.sortedCols.end(); it++)
     308            uint32_t i=0;
     309            for (auto it=fTable.sortedCols.begin(); it!=fTable.sortedCols.end(); it++, i++)
    269310            {
    270311                char *buffer = fBuffer.data() + it->offset; // pointer to column (destination buffer)
    271312
    272                 switch (it->comp)
     313                switch (fColumnOrdering[i])
    273314                {
    274                 case UNCOMPRESSED:
    275                 case SMOOTHMAN:
    276                     // regular, "semi-transposed" copy
    277                     for (char *dest=buffer; dest<buffer+thisRoundNumRows*fTable.bytes_per_row; dest+=fTable.bytes_per_row) // row-by-row
    278                     {
    279                         memcpy(dest, src, it->bytes);
    280                         src += it->bytes;  // next column
    281                     }
    282                     break;
    283 
    284                 default:
    285                     // transposed copy
    286                     for (char *elem=buffer; elem<buffer+it->bytes; elem+=it->size) // element-by-element (arrays)
    287                     {
    288                         for (char *dest=elem; dest<elem+thisRoundNumRows*fTable.bytes_per_row; dest+=fTable.bytes_per_row) // row-by-row
     315                    case FACT_ROW_MAJOR:
     316                        // regular, "semi-transposed" copy
     317                        for (char *dest=buffer; dest<buffer+thisRoundNumRows*fTable.bytes_per_row; dest+=fTable.bytes_per_row) // row-by-row
    289318                        {
    290                             memcpy(dest, src, it->size);
    291                             src += it->size; // next element
     319                            memcpy(dest, src, it->bytes);
     320                            src += it->bytes;  // next column
    292321                        }
    293                     }
     322                    break;
     323
     324                    case FACT_COL_MAJOR:
     325                        // transposed copy
     326                        for (char *elem=buffer; elem<buffer+it->bytes; elem+=it->size) // element-by-element (arrays)
     327                        {
     328                            for (char *dest=elem; dest<elem+thisRoundNumRows*fTable.bytes_per_row; dest+=fTable.bytes_per_row) // row-by-row
     329                            {
     330                                memcpy(dest, src, it->size);
     331                                src += it->size; // next element
     332                            }
     333                        }
     334                    break;
     335                    default:
     336    #ifdef __EXCEPTIONS
     337                        throw runtime_error("Unkown column ordering scheme");
     338    #else
     339                        gLog << ___err___ << "ERROR - unkown column ordering scheme" << endl;
     340                        return;
     341    #endif
    294342                    break;
    295343                };
     
    305353    uint32_t UncompressUNCOMPRESSED(char*       dest,
    306354                                    const char* src,
    307                                     uint32_t    numRows,
    308                                     uint32_t    sizeOfElems,
    309                                     uint32_t    numRowElems)
    310     {
    311         memcpy(dest, src, numRows*sizeOfElems*numRowElems);
    312         return numRows*sizeOfElems*numRowElems;
     355                                    uint32_t    numElems,
     356                                    uint32_t    sizeOfElems)
     357    {
     358        memcpy(dest, src, numElems*sizeOfElems);
     359        return numElems*sizeOfElems;
    313360    }
    314361
    315362    // Read a bunch of data compressed with the Huffman algorithm
    316     uint32_t UncompressHUFFMAN(char*       dest,
    317                                const char* src,
    318                                uint32_t ,
    319                                uint32_t    sizeOfElems,
    320                                uint32_t    numRowElems)
    321     {
    322         if (sizeOfElems < 2)
    323         {
    324             cout << "Error, Huffman only works on shorts or longer types. (here: " << sizeOfElems << "). Aborting." << endl;
    325             return -1;
    326         }
    327 
     363    uint32_t UncompressHUFFMAN16(char*       dest,
     364                                 const char* src,
     365                                 uint32_t    numChunks)
     366    {
    328367        vector<uint16_t> uncompressed;
    329368
    330369        //read compressed sizes (one per row)
    331370        const uint32_t* compressedSizes = reinterpret_cast<const uint32_t*>(src);
    332         src += sizeof(uint32_t)*numRowElems;
     371        src += sizeof(uint32_t)*numChunks;
    333372
    334373        //uncompress the rows, one by one
    335374        uint32_t sizeWritten = 0;
    336         for (uint32_t j=0;j<numRowElems;j++)
     375        for (uint32_t j=0;j<numChunks;j++)
    337376        {
    338377            Huffman::Decode(reinterpret_cast<const unsigned char*>(src), compressedSizes[j], uncompressed);
     
    347386    }
    348387
    349     //Read a bunch of data compressed with the smoothman algorithm
    350     uint32_t UncompressSMOOTHMAN(int16_t*   dest,
    351                                  const char* src,
    352                                  uint32_t numRows,
    353                                  uint32_t sizeOfElems,
    354                                  uint32_t numRowElems)
    355     {
    356         //call huffman transposed
    357         const uint32_t sizeWritten = UncompressHUFFMAN(reinterpret_cast<char*>(dest), src, numRowElems, sizeOfElems, numRows);
    358 
     388    uint32_t UnApplySMOOTHING(int16_t*   data,
     389                              uint32_t numElems)
     390    {
    359391        //un-do the integer smoothing
    360         for (uint32_t j=2;j<numRowElems*numRows;j++)
    361             dest[j] = dest[j] + (dest[j-1]+dest[j-2])/2;
    362 
    363         return sizeWritten;
    364     }
    365 
     392        for (uint32_t j=2;j<numElems;j++)
     393            data[j] = data[j] + (data[j-1]+data[j-2])/2;
     394
     395        return numElems*sizeof(uint16_t);
     396    }
    366397    // Data has been read from disk. Uncompress it !
    367398    void UncompressBuffer(const uint32_t &catalogCurrentRow, const uint32_t &thisRoundNumRows)
     
    377408
    378409            //get the compression flag
    379             const int64_t compressedOffset = fTileOffsets[catalogCurrentRow][i];//fCatalog[catalogCurrentRow][i].second - fCatalog[catalogCurrentRow][0].second;
    380             const char    compressedFlag   = fCompressedBuffer[compressedOffset];
    381 
    382             //#define COMPRESSED_FLAG 0x1
    383             //#define UNCOMPRESSED_FLAG 0x0
    384 
    385             const char *src = fCompressedBuffer.data()+compressedOffset+1;
    386 
    387             //if this bunch of data is not compressed, modify the compression flag
    388             const uint32_t compression = compressedFlag==0 ? UNCOMPRESSED : col.comp;
    389             switch (compression)
    390             {
    391                 case UNCOMPRESSED:
    392                     dest += UncompressUNCOMPRESSED(dest, src, thisRoundNumRows, col.size, col.num);
    393                     break;
    394 
    395                 case SMOOTHMAN:
    396                     dest += UncompressSMOOTHMAN(reinterpret_cast<int16_t*>(dest), src, thisRoundNumRows, col.size, col.num);
    397                     break;
    398 
    399                 default:
    400                     ;
     410            const int64_t compressedOffset = fTileOffsets[catalogCurrentRow][i];
     411
     412            BlockHeader* head = reinterpret_cast<BlockHeader*>(&fCompressedBuffer[compressedOffset]);
     413
     414            fColumnOrdering[i] = head->ordering;
     415
     416            uint32_t numRows = (head->ordering==FACT_ROW_MAJOR) ? thisRoundNumRows : col.num;
     417            uint32_t numCols = (head->ordering==FACT_COL_MAJOR) ? thisRoundNumRows : col.num;
     418
     419            const char *src = fCompressedBuffer.data()+compressedOffset+sizeof(BlockHeader)+sizeof(uint16_t)*head->numProcs;
     420
     421            for (uint32_t j=head->numProcs;j != 0; j--)
     422            {
     423                uint32_t sizeWritten=0;
     424
     425                switch (head->processings[j-1])
     426                {
     427                    case FACT_RAW:
     428                            if (head->numProcs == 1)
     429                                sizeWritten = UncompressUNCOMPRESSED(dest, src, numRows*numCols, col.size);
     430                    break;
     431                    case FACT_SMOOTHING:
     432                            sizeWritten = UnApplySMOOTHING(reinterpret_cast<int16_t*>(dest), numRows*numCols);
     433                    break;
     434                    case FACT_HUFFMAN16:
     435                            sizeWritten = UncompressHUFFMAN16(dest, src, numRows);
     436                    break;
     437                    default:
     438#ifdef __EXCEPTIONS
     439                    throw runtime_error("Unknown processing applied to data. Aborting");
     440#else
     441                    gLog << ___err___ << "ERROR - Unknown processing applied to data. Aborting" << endl;
     442                    return;
     443#endif
     444                    break;
     445                }
     446                //increment destination counter only when processing done.
     447                if (j==1) dest+= sizeWritten;
    401448            }
    402449        }
Note: See TracChangeset for help on using the changeset viewer.