Index: trunk/Mars/mbase/MFits.h
===================================================================
--- trunk/Mars/mbase/MFits.h	(revision 11549)
+++ 	(revision )
@@ -1,562 +1,0 @@
-#ifndef MARS_MFits
-#define MARS_MFits
-
-#include <map>
-#include <string>
-#include <sstream>
-#include <algorithm>
-
-#include <stdint.h>
-
-#ifdef __CINT__
-#define off_t size_t
-#endif
-
-#ifndef __MARS__
-#include <vector>
-#include <fstream>
-#include <iomanip>
-#include <iostream>
-#define MZlib ifstream
-#define gLog cerr
-#define err ""
-#define all ""
-#else
-#include "MZlib.h"
-
-#include "MLog.h"
-#include "MLogManip.h"
-#endif
-
-using namespace std;
-
-class MFits : public MZlib
-{
-public:
-    struct Entry
-    {
-        char   type;
-        string value;
-        string comment;
-
-        template<typename T>
-            T Get() const
-        {
-            T t;
-
-            istringstream str(value);
-            str >> t;
-
-            return t;
-        }
-    };
-
-    struct Table
-    {
-        off_t offset;
-
-        string name;
-        size_t bytes_per_row;
-        size_t num_rows;
-        size_t num_cols;
-
-        struct Column
-        {
-            size_t offset;
-            size_t num;
-            size_t size;
-            char   type;
-        };
-
-        typedef map<string, Entry>  Keys;
-        typedef map<string, Column> Columns;
-
-        Columns cols;
-        Keys    keys;
-
-        string Trim(const string &str, char c=' ') const
-        {
-            // Trim Both leading and trailing spaces
-            const size_t start = str.find_first_not_of(c); // Find the first character position after excluding leading blank spaces
-            const size_t end   = str.find_last_not_of(c);  // Find the first character position from reverse af
-
-            // if all spaces or empty return an empty string
-            if (string::npos==start || string::npos==end)
-                return string();
-
-            return str.substr(start, end-start+1);
-        }
-
-        bool Check(const string &key, char type, const string &value="") const
-        {
-            const Keys::const_iterator it = keys.find(key);
-            if (it==keys.end())
-            {
-                gLog << err << "ERROR - Key '" << key << "' not found." << endl;
-                return false;
-            }
-
-            if (it->second.type!=type)
-            {
-                gLog << err << "ERROR - Wrong type for key '" << key << "': expected " << type << ", found " << it->second.type << "." << endl;
-                return false;
-            }
-
-            if (!value.empty() && it->second.value!=value)
-            {
-                gLog << err << "ERROR - Wrong value for key '" << key << "': expected " << value << ", found " << it->second.value << "." << endl;
-                return false;
-            }
-
-            return true;
-        }
-
-        Keys ParseBlock(const vector<string> &vec) const
-        {
-            map<string,Entry> rc;
-
-            for (unsigned int i=0; i<vec.size(); i++)
-            {
-                const string key = Trim(vec[i].substr(0,8));
-                // Keywords without a value, like COMMENT / HISTORY
-                if (vec[i].substr(8,2)!="= ")
-                    continue;
-
-                char type = 0;
-
-                string com;
-                string val = Trim(vec[i].substr(10));
-                if (val[0]=='\'')
-                {
-                    // First skip all '' in the string
-                    size_t p = 1;
-                    while (1)
-                    {
-                        const size_t pp = val.find_first_of('\'', p);
-                        if (pp==string::npos)
-                            break;
-
-                        p = val[pp+1]=='\'' ? pp+2 : pp+1;
-                    }
-
-                    // Now find the comment
-                    const size_t ppp = val.find_first_of('/', p);
-
-                    // Set value, comment and type
-                    com  = ppp==string::npos ? "" : Trim(val.substr(ppp+1));
-                    val  = Trim(val.substr(1, p-2));
-                    type = 'T';
-                }
-                else
-                {
-                    const size_t p = val.find_first_of('/');
-
-                    com = Trim(val.substr(p+2));
-                    val = Trim(val.substr(0, p));
-
-                    if (val.empty() || val.find_first_of('T')!=string::npos || val.find_first_of('F')!=string::npos)
-                        type = 'B';
-                    else
-                        type = val.find_last_of('.')==string::npos ? 'I' : 'F';
-                }
-
-                const Entry e = { type, val, com };
-                rc[key] = e;
-            }
-
-            return rc;
-        }
-
-        Table() : offset(0) { }
-        Table(const vector<string> &vec, off_t off) :
-            offset(off), keys(ParseBlock(vec))
-        {
-            if (!Check("XTENSION", 'T', "BINTABLE") ||
-                !Check("NAXIS",    'I', "2")        ||
-                !Check("BITPIX",   'I', "8")        ||
-                !Check("PCOUNT",   'I', "0")        ||
-                !Check("GCOUNT",   'I', "1")        ||
-                !Check("EXTNAME",  'T')             ||
-                !Check("NAXIS1",   'I')             ||
-                !Check("NAXIS2",   'I')             ||
-                !Check("TFIELDS",  'I'))
-                return;
-
-            bytes_per_row = Get<size_t>("NAXIS1");
-            num_rows      = Get<size_t>("NAXIS2");
-            num_cols      = Get<size_t>("TFIELDS");
-
-            size_t bytes = 0;
-            for (size_t i=1; i<=num_cols; i++)
-            {
-                ostringstream num;
-                num << i;
-
-                if (!Check("TTYPE"+num.str(), 'T') ||
-                    !Check("TFORM"+num.str(), 'T'))
-                    return;
-
-                const string name = Get<string>("TTYPE"+num.str());
-                const string fmt  = Get<string>("TFORM"+num.str());
-
-                istringstream in(fmt);
-                int n = 0;
-                in >> n;
-                if (!in)
-                    n = 1;
-
-                const char type = fmt[fmt.length()-1];
-
-                size_t size = 0;
-                switch (type)
-                {
-                    // We could use negative values to mark floats
-                    // otheriwse we could just cast them to int64_t?
-                case 'L': size = 1; break; // logical
-                // case 'X': size = n; break; // bits (n=number of bytes needed to contain all bits)
-                case 'B': size = 1; break; // byte
-                case 'I': size = 2; break; // short
-                case 'J': size = 4; break; // int
-                case 'K': size = 8; break; // long long
-                case 'E': size = 4; break; // float
-                case 'D': size = 8; break; // double
-                // case 'C': size =  8; break; // complex float
-                // case 'M': size = 16; break; // complex double
-                // case 'P': size =  8; break; // array descriptor (32bit)
-                // case 'Q': size = 16; break; // array descriptor (64bit)
-                default:
-                    gLog << err << "FITS format '" << fmt << "' not yet supported." << endl;
-                    continue;
-                }
-
-                const Table::Column col = { bytes, n, size, type };
-
-                cols[name] = col;
-                bytes  += n*size;
-            }
-
-            if (bytes!=bytes_per_row)
-            {
-                gLog << err << "Column size mismatch" << endl;
-                return;
-            }
-
-            name = Get<string>("EXTNAME");
-        }
-
-        void PrintKeys() const
-        {
-            for (Keys::const_iterator it=keys.begin(); it!=keys.end(); it++)
-                gLog << all << setw(2) << it->second.type << '|' << it->first << '=' << it->second.value << '/' << it->second.comment << '|' << endl;
-        }
-
-        void PrintColumns() const
-        {
-            typedef map<pair<size_t, string>, Column> Sorted;
-
-            Sorted sorted;
-
-            for (Columns::const_iterator it=cols.begin(); it!=cols.end(); it++)
-                sorted[make_pair(it->second.offset, it->first)] = it->second;
-
-            for (Sorted::const_iterator it=sorted.begin(); it!=sorted.end(); it++)
-                gLog << all << setw(6) << it->second.offset << "| " << it->first.second << " [" << it->second.type << '=' << it->second.num << 'x' << it->second.size << "byte]" << endl;
-        }
-
-        operator bool() const { return !name.empty(); }
-
-        bool HasKey(const string &key) const
-        {
-            return keys.find(key)!=keys.end();
-        }
-
-            // Values of keys are always signed
-        template<typename T>
-            T Get(const string &key) const
-        {
-            const ::map<string,Entry>::const_iterator it = keys.find(key);
-            if (it==keys.end())
-            {
-                gLog << err << "Key '" << key << "' not found." << endl;
-                return 0;
-            }
-            return it->second.Get<T>();
-        }
-
-        size_t GetN(const string &key) const
-        {
-            const Columns::const_iterator it = cols.find(key);
-            return it==cols.end() ? 0 : it->second.num;
-        }
-    };
-
-private:
-    Table fTable;
-
-    typedef pair<void*, Table::Column> Address;
-    typedef vector<Address> Addresses;
-    //map<void*, Table::Column> fAddresses;
-    Addresses fAddresses;
-
-    vector<char> fBufferRow;
-    vector<char> fBufferDat;
-
-    size_t fRow;
-
-    vector<string> ReadBlock(vector<string> &vec)
-    {
-        bool end = false;
-        for (int i=0; i<36; i++)
-        {
-            char c[81];
-            c[80] = 0;
-            read(c, 80);
-            if (!good())
-                break;
-
-            if (c[0]==0)
-                return vector<string>();
-
-            string str(c);
-
-            //if (!str.empty())
-            //    cout << setw(2) << i << "|" << str << "|" << endl;
-
-            if (str=="END                                                                             ")
-                end = true;
-
-            if (end)
-                continue;
-
-            vec.push_back(str);
-        }
-
-        return vec;
-    }
-
-    string Compile(const string &key, int16_t i=-1)
-    {
-        if (i<0)
-            return key;
-
-        ostringstream str;
-        str << key << i;
-        return str.str();
-    }
-
-public:
-    MFits(const string &fname) : MZlib(fname.c_str())
-    {
-        while (good())
-        {
-            vector<string> block;
-            while (1)
-            {
-                // FIXME: Set limit on memory consumption
-                ReadBlock(block);
-                if (!good())
-                {
-                    gLog << err << "FITS file corrupted." << endl;
-                    clear(rdstate()|ios::badbit);
-                    return;
-                }
-
-                if (block.size()%36)
-                    break;
-            }
-
-            if (block.size()==0)
-                break;
-
-            if (block[0].substr(0, 9)=="SIMPLE  =")
-                continue;
-
-            if (block[0].substr(0, 9)=="XTENSION=")
-            {
-                // FIXME: Check for table name
-
-                fTable = Table(block, tellg());
-                fRow   = (size_t)-1;
-
-                fTable.PrintKeys();
-
-                if (!fTable)
-                    return;
-
-                if (!fTable.Check("TELESCOP", 'T', "FACT"))
-                {
-                    gLog << err << "Not a valid FACT fits file." << endl;
-                    return;
-                }
-
-                fBufferRow.resize(fTable.bytes_per_row);
-                fBufferDat.resize(fTable.bytes_per_row);
-
-                /*
-                 // Next table should start at:
-                 const size_t size = fTable.bytes_per_row*fTable.num_rows;
-                 const size_t blks = size/(36*80);
-                 const size_t rest = size%(36*80);
-
-                seekg((blks+(rest>0?1:0))*(36*80), ios::cur);
-                if (!good())
-                     gLog << err << "File seems to be incomplete (less data than expected from header)." << endl;
-
-                fRow   = fTable.num_rows;
-                 */
-
-                break;
-            }
-        }
-    }
-
-    void ReadRow(size_t row)
-    {
-        // if (row!=fRow+1) // Fast seeking is ensured by MZlib
-        seekg(fTable.offset+row*fTable.bytes_per_row);
-
-        fRow = row;
-
-        read(fBufferRow.data(), fBufferRow.size());
-        //fin.clear(fin.rdstate()&~ios::eofbit);
-    }
-
-    template<size_t N>
-        void revcpy(char *dest, const char *src, int num)
-    {
-        const char *end = src + num*N;
-        for (const char *ptr = src; ptr<end; ptr+=N, dest+=N)
-            reverse_copy(ptr, ptr+N, dest);
-    }
-
-    bool GetRow(size_t row)
-    {
-        ReadRow(row);
-        if (!good())
-            return good();
-
-        for (Addresses::const_iterator it=fAddresses.begin(); it!=fAddresses.end(); it++)
-        {
-            const Table::Column &c = it->second;
-
-            const char *src = fBufferRow.data() + c.offset;
-            char *dest = reinterpret_cast<char*>(it->first);
-
-            // Let the compiler do some optimization by
-            // knowing the 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;
-            }
-        }
-
-        return good();
-    }
-
-    bool GetNextRow()
-    {
-        return GetRow(fRow+1);
-    }
-
-    bool SkipNextRow()
-    {
-        seekg(fTable.offset+(++fRow)*fTable.bytes_per_row);
-        return good();
-    }
-
-    static bool Compare(const Address &p1, const Address &p2)
-    {
-        return p1.first>p2.first;
-    }
-
-    template<typename T>
-    bool SetPtrAddress(const string &name, T *ptr, size_t cnt)
-    {
-        if (fTable.cols.count(name)==0)
-        {
-            gLog << err <<"SetPtrAddress('" << name << "') - Column not found." << endl;
-            return false;
-        }
-
-        if (sizeof(T)!=fTable.cols[name].size)
-        {
-            gLog << err << "SetPtrAddress('" << name << "') - Element size mismatch: expected "
-                << fTable.cols[name].size << " from header, got " << sizeof(T) << endl;
-
-            return false;
-        }
-
-        if (cnt!=fTable.cols[name].num)
-        {
-            gLog << err << "SetPtrAddress('" << name << "') - Element count mismatch: expected "
-                << fTable.cols[name].num << " from header, got " << cnt << endl;
-
-            return false;
-        }
-
-        // if (fAddresses.count(ptr)>0)
-        //     gLog << warn << "SetPtrAddress('" << name << "') - Pointer " << ptr << " already assigned." << endl;
-
-        //fAddresses[ptr] = fTable.cols[name];
-        fAddresses.push_back(make_pair(ptr, fTable.cols[name]));
-        sort(fAddresses.begin(), fAddresses.end(), Compare);
-        return true;
-    }
-
-    template<class T>
-    bool SetRefAddress(const string &name, T &ptr)
-    {
-        return SetPtrAddress(name, &ptr, sizeof(ptr)/sizeof(T));
-    }
-
-    template<typename T>
-    bool SetVecAddress(const string &name, vector<T> &vec)
-    {
-        return SetPtrAddress(name, vec.data(), vec.size());
-    }
-
-    template<typename T>
-        T Get(const string &key) const
-    {
-        return fTable.Get<T>(key);
-    }
-
-    bool SetPtrAddress(const string &name, void *ptr)
-    {
-        if (fTable.cols.count(name)==0)
-        {
-            gLog << err <<"SetPtrAddress('" << name << "') - Column not found." << endl;
-            return false;
-        }
-
-        // if (fAddresses.count(ptr)>0)
-        //     gLog << warn << "SetPtrAddress('" << name << "') - Pointer " << ptr << " already assigned." << endl;
-
-        //fAddresses[ptr] = fTable.cols[name];
-        fAddresses.push_back(make_pair(ptr, fTable.cols[name]));
-        sort(fAddresses.begin(), fAddresses.end(), Compare);
-        return true;
-    }
-
-    bool     HasKey(const string &key) const { return fTable.HasKey(key); }
-    int64_t  GetInt(const string &key) const { return fTable.Get<int64_t>(key); }
-    uint64_t GetUInt(const string &key) const { return fTable.Get<uint64_t>(key); }
-    double   GetFloat(const string &key) const { return fTable.Get<double>(key); }
-
-    size_t GetN(const string &key) const
-    {
-        return fTable.GetN(key);
-    }
-
-    size_t GetNumRows() const { return fTable.num_rows; }
-    size_t GetRow() const { return fRow; }
-
-    operator bool() const { return fTable && fTable.offset!=0; }
-
-    void PrintKeys() const { fTable.PrintKeys(); }
-    void PrintColumns() const { fTable.PrintColumns(); }
-};
-#endif
Index: trunk/Mars/mbase/fits.h
===================================================================
--- trunk/Mars/mbase/fits.h	(revision 11550)
+++ trunk/Mars/mbase/fits.h	(revision 11550)
@@ -0,0 +1,576 @@
+#ifndef MARS_MFits
+#define MARS_MFits
+
+#include <map>
+#include <string>
+#include <sstream>
+#include <algorithm>
+
+#include <stdint.h>
+
+#ifdef __CINT__
+#define off_t size_t
+#endif
+
+#ifndef __MARS__
+#include <vector>
+#include <fstream>
+#include <iomanip>
+#include <iostream>
+#define MZlib ifstream
+#define gLog cerr
+#define err ""
+#define all ""
+#else
+#include "MZlib.h"
+
+#include "MLog.h"
+#include "MLogManip.h"
+#endif
+
+using namespace std;
+
+class MFits : public MZlib
+{
+public:
+    struct Entry
+    {
+        char   type;
+        string value;
+        string comment;
+
+        template<typename T>
+            T Get() const
+        {
+            T t;
+
+            istringstream str(value);
+            str >> t;
+
+            return t;
+        }
+    };
+
+    struct Table
+    {
+        off_t offset;
+
+        string name;
+        size_t bytes_per_row;
+        size_t num_rows;
+        size_t num_cols;
+
+        struct Column
+        {
+            size_t offset;
+            size_t num;
+            size_t size;
+            char   type;
+        };
+
+        typedef map<string, Entry>  Keys;
+        typedef map<string, Column> Columns;
+
+        Columns cols;
+        Keys    keys;
+
+        string Trim(const string &str, char c=' ') const
+        {
+            // Trim Both leading and trailing spaces
+            const size_t start = str.find_first_not_of(c); // Find the first character position after excluding leading blank spaces
+            const size_t end   = str.find_last_not_of(c);  // Find the first character position from reverse af
+
+            // if all spaces or empty return an empty string
+            if (string::npos==start || string::npos==end)
+                return string();
+
+            return str.substr(start, end-start+1);
+        }
+
+        bool Check(const string &key, char type, const string &value="") const
+        {
+            const Keys::const_iterator it = keys.find(key);
+            if (it==keys.end())
+            {
+                gLog << err << "ERROR - Key '" << key << "' not found." << endl;
+                return false;
+            }
+
+            if (it->second.type!=type)
+            {
+                gLog << err << "ERROR - Wrong type for key '" << key << "': expected " << type << ", found " << it->second.type << "." << endl;
+                return false;
+            }
+
+            if (!value.empty() && it->second.value!=value)
+            {
+                gLog << err << "ERROR - Wrong value for key '" << key << "': expected " << value << ", found " << it->second.value << "." << endl;
+                return false;
+            }
+
+            return true;
+        }
+
+        Keys ParseBlock(const vector<string> &vec) const
+        {
+            map<string,Entry> rc;
+
+            for (unsigned int i=0; i<vec.size(); i++)
+            {
+                const string key = Trim(vec[i].substr(0,8));
+                // Keywords without a value, like COMMENT / HISTORY
+                if (vec[i].substr(8,2)!="= ")
+                    continue;
+
+                char type = 0;
+
+                string com;
+                string val = Trim(vec[i].substr(10));
+                if (val[0]=='\'')
+                {
+                    // First skip all '' in the string
+                    size_t p = 1;
+                    while (1)
+                    {
+                        const size_t pp = val.find_first_of('\'', p);
+                        if (pp==string::npos)
+                            break;
+
+                        p = val[pp+1]=='\'' ? pp+2 : pp+1;
+                    }
+
+                    // Now find the comment
+                    const size_t ppp = val.find_first_of('/', p);
+
+                    // Set value, comment and type
+                    com  = ppp==string::npos ? "" : Trim(val.substr(ppp+1));
+                    val  = Trim(val.substr(1, p-2));
+                    type = 'T';
+                }
+                else
+                {
+                    const size_t p = val.find_first_of('/');
+
+                    com = Trim(val.substr(p+2));
+                    val = Trim(val.substr(0, p));
+
+                    if (val.empty() || val.find_first_of('T')!=string::npos || val.find_first_of('F')!=string::npos)
+                        type = 'B';
+                    else
+                        type = val.find_last_of('.')==string::npos ? 'I' : 'F';
+                }
+
+                const Entry e = { type, val, com };
+                rc[key] = e;
+            }
+
+            return rc;
+        }
+
+        Table() : offset(0) { }
+        Table(const vector<string> &vec, off_t off) :
+            offset(off), keys(ParseBlock(vec))
+        {
+            if (!Check("XTENSION", 'T', "BINTABLE") ||
+                !Check("NAXIS",    'I', "2")        ||
+                !Check("BITPIX",   'I', "8")        ||
+                !Check("PCOUNT",   'I', "0")        ||
+                !Check("GCOUNT",   'I', "1")        ||
+                !Check("EXTNAME",  'T')             ||
+                !Check("NAXIS1",   'I')             ||
+                !Check("NAXIS2",   'I')             ||
+                !Check("TFIELDS",  'I'))
+                return;
+
+            bytes_per_row = Get<size_t>("NAXIS1");
+            num_rows      = Get<size_t>("NAXIS2");
+            num_cols      = Get<size_t>("TFIELDS");
+
+            size_t bytes = 0;
+            for (size_t i=1; i<=num_cols; i++)
+            {
+                ostringstream num;
+                num << i;
+
+                if (!Check("TTYPE"+num.str(), 'T') ||
+                    !Check("TFORM"+num.str(), 'T'))
+                    return;
+
+                const string name = Get<string>("TTYPE"+num.str());
+                const string fmt  = Get<string>("TFORM"+num.str());
+
+                istringstream in(fmt);
+                int n = 0;
+                in >> n;
+                if (!in)
+                    n = 1;
+
+                const char type = fmt[fmt.length()-1];
+
+                size_t size = 0;
+                switch (type)
+                {
+                    // We could use negative values to mark floats
+                    // otheriwse we could just cast them to int64_t?
+                case 'L': size = 1; break; // logical
+                // case 'X': size = n; break; // bits (n=number of bytes needed to contain all bits)
+                case 'B': size = 1; break; // byte
+                case 'I': size = 2; break; // short
+                case 'J': size = 4; break; // int
+                case 'K': size = 8; break; // long long
+                case 'E': size = 4; break; // float
+                case 'D': size = 8; break; // double
+                // case 'C': size =  8; break; // complex float
+                // case 'M': size = 16; break; // complex double
+                // case 'P': size =  8; break; // array descriptor (32bit)
+                // case 'Q': size = 16; break; // array descriptor (64bit)
+                default:
+                    gLog << err << "FITS format '" << fmt << "' not yet supported." << endl;
+                    continue;
+                }
+
+                const Table::Column col = { bytes, n, size, type };
+
+                cols[name] = col;
+                bytes  += n*size;
+            }
+
+            if (bytes!=bytes_per_row)
+            {
+                gLog << err << "Column size mismatch" << endl;
+                return;
+            }
+
+            name = Get<string>("EXTNAME");
+        }
+
+        void PrintKeys() const
+        {
+            for (Keys::const_iterator it=keys.begin(); it!=keys.end(); it++)
+                gLog << all << setw(2) << it->second.type << '|' << it->first << '=' << it->second.value << '/' << it->second.comment << '|' << endl;
+        }
+
+        void PrintColumns() const
+        {
+            typedef map<pair<size_t, string>, Column> Sorted;
+
+            Sorted sorted;
+
+            for (Columns::const_iterator it=cols.begin(); it!=cols.end(); it++)
+                sorted[make_pair(it->second.offset, it->first)] = it->second;
+
+            for (Sorted::const_iterator it=sorted.begin(); it!=sorted.end(); it++)
+                gLog << all << setw(6) << it->second.offset << "| " << it->first.second << " [" << it->second.type << '=' << it->second.num << 'x' << it->second.size << "byte]" << endl;
+        }
+
+        operator bool() const { return !name.empty(); }
+
+        bool HasKey(const string &key) const
+        {
+            return keys.find(key)!=keys.end();
+        }
+
+            // Values of keys are always signed
+        template<typename T>
+            T Get(const string &key) const
+        {
+            const ::map<string,Entry>::const_iterator it = keys.find(key);
+            if (it==keys.end())
+            {
+                gLog << err << "Key '" << key << "' not found." << endl;
+                return 0;
+            }
+            return it->second.Get<T>();
+        }
+
+        size_t GetN(const string &key) const
+        {
+            const Columns::const_iterator it = cols.find(key);
+            return it==cols.end() ? 0 : it->second.num;
+        }
+    };
+
+private:
+    Table fTable;
+
+    typedef pair<void*, Table::Column> Address;
+    typedef vector<Address> Addresses;
+    //map<void*, Table::Column> fAddresses;
+    Addresses fAddresses;
+
+    vector<char> fBufferRow;
+    vector<char> fBufferDat;
+
+    size_t fRow;
+
+    vector<string> ReadBlock(vector<string> &vec)
+    {
+        bool end = false;
+        for (int i=0; i<36; i++)
+        {
+            char c[81];
+            c[80] = 0;
+            read(c, 80);
+            if (!good())
+                break;
+
+            if (c[0]==0)
+                return vector<string>();
+
+            string str(c);
+
+            //if (!str.empty())
+            //    cout << setw(2) << i << "|" << str << "|" << endl;
+
+            if (str=="END                                                                             ")
+                end = true;
+
+            if (end)
+                continue;
+
+            vec.push_back(str);
+        }
+
+        return vec;
+    }
+
+    string Compile(const string &key, int16_t i=-1)
+    {
+        if (i<0)
+            return key;
+
+        ostringstream str;
+        str << key << i;
+        return str.str();
+    }
+
+public:
+    MFits(const string &fname) : MZlib(fname.c_str())
+    {
+        char simple[6];
+        read(simple, 6);
+        if (!good())
+            return;
+
+        if (memcmp(simple, "SIMPLE", 6))
+        {
+            gLog << err << "File is not a FITS file." << endl;
+            clear(rdstate()|ios::badbit);
+            return;
+        }
+
+        seekpos(0);
+
+        while (good())
+        {
+            vector<string> block;
+            while (1)
+            {
+                // FIXME: Set limit on memory consumption
+                ReadBlock(block);
+                if (!good())
+                {
+                    gLog << err << "FITS file corrupted." << endl;
+                    clear(rdstate()|ios::badbit);
+                    return;
+                }
+
+                if (block.size()%36)
+                    break;
+            }
+
+            if (block.size()==0)
+                break;
+
+            if (block[0].substr(0, 9)=="SIMPLE  =")
+                continue;
+
+            if (block[0].substr(0, 9)=="XTENSION=")
+            {
+                // FIXME: Check for table name
+
+                fTable = Table(block, tellg());
+                fRow   = (size_t)-1;
+
+                fTable.PrintKeys();
+
+                if (!fTable)
+                    return;
+
+                if (!fTable.Check("TELESCOP", 'T', "FACT"))
+                {
+                    gLog << err << "Not a valid FACT fits file." << endl;
+                    return;
+                }
+
+                fBufferRow.resize(fTable.bytes_per_row);
+                fBufferDat.resize(fTable.bytes_per_row);
+
+                /*
+                 // Next table should start at:
+                 const size_t size = fTable.bytes_per_row*fTable.num_rows;
+                 const size_t blks = size/(36*80);
+                 const size_t rest = size%(36*80);
+
+                seekg((blks+(rest>0?1:0))*(36*80), ios::cur);
+                if (!good())
+                     gLog << err << "File seems to be incomplete (less data than expected from header)." << endl;
+
+                fRow   = fTable.num_rows;
+                 */
+
+                break;
+            }
+        }
+    }
+
+    void ReadRow(size_t row)
+    {
+        // if (row!=fRow+1) // Fast seeking is ensured by MZlib
+        seekg(fTable.offset+row*fTable.bytes_per_row);
+
+        fRow = row;
+
+        read(fBufferRow.data(), fBufferRow.size());
+        //fin.clear(fin.rdstate()&~ios::eofbit);
+    }
+
+    template<size_t N>
+        void revcpy(char *dest, const char *src, int num)
+    {
+        const char *end = src + num*N;
+        for (const char *ptr = src; ptr<end; ptr+=N, dest+=N)
+            reverse_copy(ptr, ptr+N, dest);
+    }
+
+    bool GetRow(size_t row)
+    {
+        ReadRow(row);
+        if (!good())
+            return good();
+
+        for (Addresses::const_iterator it=fAddresses.begin(); it!=fAddresses.end(); it++)
+        {
+            const Table::Column &c = it->second;
+
+            const char *src = fBufferRow.data() + c.offset;
+            char *dest = reinterpret_cast<char*>(it->first);
+
+            // Let the compiler do some optimization by
+            // knowing the 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;
+            }
+        }
+
+        return good();
+    }
+
+    bool GetNextRow()
+    {
+        return GetRow(fRow+1);
+    }
+
+    bool SkipNextRow()
+    {
+        seekg(fTable.offset+(++fRow)*fTable.bytes_per_row);
+        return good();
+    }
+
+    static bool Compare(const Address &p1, const Address &p2)
+    {
+        return p1.first>p2.first;
+    }
+
+    template<typename T>
+    bool SetPtrAddress(const string &name, T *ptr, size_t cnt)
+    {
+        if (fTable.cols.count(name)==0)
+        {
+            gLog << err <<"SetPtrAddress('" << name << "') - Column not found." << endl;
+            return false;
+        }
+
+        if (sizeof(T)!=fTable.cols[name].size)
+        {
+            gLog << err << "SetPtrAddress('" << name << "') - Element size mismatch: expected "
+                << fTable.cols[name].size << " from header, got " << sizeof(T) << endl;
+
+            return false;
+        }
+
+        if (cnt!=fTable.cols[name].num)
+        {
+            gLog << err << "SetPtrAddress('" << name << "') - Element count mismatch: expected "
+                << fTable.cols[name].num << " from header, got " << cnt << endl;
+
+            return false;
+        }
+
+        // if (fAddresses.count(ptr)>0)
+        //     gLog << warn << "SetPtrAddress('" << name << "') - Pointer " << ptr << " already assigned." << endl;
+
+        //fAddresses[ptr] = fTable.cols[name];
+        fAddresses.push_back(make_pair(ptr, fTable.cols[name]));
+        sort(fAddresses.begin(), fAddresses.end(), Compare);
+        return true;
+    }
+
+    template<class T>
+    bool SetRefAddress(const string &name, T &ptr)
+    {
+        return SetPtrAddress(name, &ptr, sizeof(ptr)/sizeof(T));
+    }
+
+    template<typename T>
+    bool SetVecAddress(const string &name, vector<T> &vec)
+    {
+        return SetPtrAddress(name, vec.data(), vec.size());
+    }
+
+    template<typename T>
+        T Get(const string &key) const
+    {
+        return fTable.Get<T>(key);
+    }
+
+    bool SetPtrAddress(const string &name, void *ptr)
+    {
+        if (fTable.cols.count(name)==0)
+        {
+            gLog << err <<"SetPtrAddress('" << name << "') - Column not found." << endl;
+            return false;
+        }
+
+        // if (fAddresses.count(ptr)>0)
+        //     gLog << warn << "SetPtrAddress('" << name << "') - Pointer " << ptr << " already assigned." << endl;
+
+        //fAddresses[ptr] = fTable.cols[name];
+        fAddresses.push_back(make_pair(ptr, fTable.cols[name]));
+        sort(fAddresses.begin(), fAddresses.end(), Compare);
+        return true;
+    }
+
+    bool     HasKey(const string &key) const { return fTable.HasKey(key); }
+    int64_t  GetInt(const string &key) const { return fTable.Get<int64_t>(key); }
+    uint64_t GetUInt(const string &key) const { return fTable.Get<uint64_t>(key); }
+    double   GetFloat(const string &key) const { return fTable.Get<double>(key); }
+
+    size_t GetN(const string &key) const
+    {
+        return fTable.GetN(key);
+    }
+
+    size_t GetNumRows() const { return fTable.num_rows; }
+    size_t GetRow() const { return fRow; }
+
+    operator bool() const { return fTable && fTable.offset!=0; }
+
+    void PrintKeys() const { fTable.PrintKeys(); }
+    void PrintColumns() const { fTable.PrintColumns(); }
+};
+#endif