source: trunk/FACT++/src/fitsdump.cc @ 10837

Last change on this file since 10837 was 10837, checked in by tbretz, 8 years ago
Some performance improvements; added some comments with code to fill just a histogram or a graph - maybe for later use.
File size: 20.5 KB
Line 
1//****************************************************************
2/** @class FitsDumper
3
4  @brief Dumps contents of fits tables to stdout or a file
5
6 */
7 //****************************************************************
8#include "Configuration.h"
9
10#include <map>
11#include <fstream>
12
13#include <CCfits/CCfits>
14
15using namespace std;
16
17class FitsDumper
18{
19public:
20    FitsDumper();
21    ~FitsDumper();
22
23private:
24   
25    CCfits::FITS  *fFile;                 /// FITS pointer
26    CCfits::Table *fTable;                /// Table pointer
27    map<string, CCfits::Column*> fColMap; /// map between the column names and their CCfits objects
28
29    // Convert CCfits::ValueType into a human readable string
30    string ValueTypeToStr(CCfits::ValueType type) const;
31
32    // Convert CCfits::ValueType into a number of associated bytes
33    int    ValueTypeToSize(CCfits::ValueType type) const;
34
35    /// Calculate the buffer size required to read a row of the fits table, as well as the offsets to each column
36    vector<int> CalculateOffsets() const;
37
38    template<class T>
39        T PtrToValue(const unsigned char* &ptr) const;
40//    template<class T>
41//        double PtrToDouble(const unsigned char *ptr) const;
42//    double PtrToDouble(const unsigned char *ptr, CCfits::ValueType type) const;
43
44    /// Write a single row of the selected data
45    int  WriteRow(ostream &, const vector<CCfits::Column*> &, const vector<int> &, unsigned char *) const;
46
47    bool OpenFile(const string &);        /// Open a file
48    bool OpenTable(const string &);       /// Open a table
49
50    /// Lists all columns of an open file
51    void List();                         
52    void ListHeader();
53    void ListKeywords(ostream &);
54
55    /// Perform the dumping, based on the current dump list
56    bool Dump(const string &, const vector<string> &list, int);
57
58//    bool Plot(const vector<string> &list);
59
60public:
61    ///Configures the fitsLoader from the config file and/or command arguments.
62    int ExecConfig(Configuration& conf);
63};
64
65// --------------------------------------------------------------------------
66//
67//! Constructor
68//! @param out
69//!        the ostream where to redirect the outputs
70//
71FitsDumper::FitsDumper() : fFile(0), fTable(0)
72{
73}
74
75// --------------------------------------------------------------------------
76//
77//! Destructor
78//
79FitsDumper::~FitsDumper()
80{
81    if (fFile)
82        delete fFile;
83}
84
85
86string FitsDumper::ValueTypeToStr(CCfits::ValueType type) const
87{
88    switch (type)
89    {
90    case CCfits::Tbyte:     return "uint8_t";
91    case CCfits::Tushort:   return "uint16_t";
92    case CCfits::Tshort:    return "int16_t";
93    case CCfits::Tuint:     return "uint32_t";
94    case CCfits::Tint:      return "int32_t";
95    case CCfits::Tulong:    return "uint32_t";
96    case CCfits::Tlong:     return "int32_t";
97    case CCfits::Tlonglong: return "int64_t";
98    case CCfits::Tfloat:    return "float";
99    case CCfits::Tdouble:   return "double";
100
101    default:
102        return "unknwown";
103    }
104}
105
106int FitsDumper::ValueTypeToSize(CCfits::ValueType type) const
107{
108    switch (type)
109    {
110    case CCfits::Tbyte:     return sizeof(uint8_t);
111    case CCfits::Tushort:
112    case CCfits::Tshort:    return sizeof(uint16_t);
113    case CCfits::Tuint:
114    case CCfits::Tint:
115    case CCfits::Tulong:
116    case CCfits::Tlong:     return sizeof(uint32_t);
117    case CCfits::Tlonglong: return sizeof(uint64_t);
118    case CCfits::Tfloat:    return sizeof(float);
119    case CCfits::Tdouble:   return sizeof(double);
120
121    default:
122        return 0;
123    }
124}
125
126template<class T>
127T FitsDumper::PtrToValue(const unsigned char* &ptr) const
128{
129    T t;
130    reverse_copy(ptr, ptr+sizeof(T), reinterpret_cast<unsigned char*>(&t));
131    ptr += sizeof(T);
132
133    return t;
134}
135/*
136template<class T>
137double FitsDumper::PtrToDouble(const unsigned char *ptr) const
138{
139    T t;
140    reverse_copy(ptr, ptr+sizeof(T), reinterpret_cast<unsigned char*>(&t));
141    return t;
142}
143
144double FitsDumper::PtrToDouble(const unsigned char *ptr, CCfits::ValueType type) const
145{
146    switch (type)
147    {
148    case CCfits::Tbyte:     return PtrToDouble<uint8_t> (ptr);
149    case CCfits::Tushort:   return PtrToDouble<uint16_t>(ptr);
150    case CCfits::Tuint:
151    case CCfits::Tulong:    return PtrToDouble<uint32_t>(ptr);
152    case CCfits::Tshort:    return PtrToDouble<int16_t> (ptr);
153    case CCfits::Tint:
154    case CCfits::Tlong:     return PtrToDouble<int32_t> (ptr);
155    case CCfits::Tlonglong: return PtrToDouble<int64_t> (ptr);
156    case CCfits::Tfloat:    return PtrToDouble<float>   (ptr);
157    case CCfits::Tdouble:   return PtrToDouble<double>  (ptr);
158
159    default:
160        throw runtime_error("Data type not implemented yet.");
161    }
162}
163*/
164
165// --------------------------------------------------------------------------
166//
167//! Writes a single row of the selected FITS data to the output file.
168//!
169//! Data type \b not yet implemented:
170//!   CCfits::Tnull, CCfits::Tbit, CCfits::Tlogical, CCfits::Tstring,
171//!   CCfits::Tcomplex, CCfits::Tdblcomplex, CCfits::VTbit,
172//!   CCfits::VTbyte, CCfits::VTlogical, CCfits::VTushort,
173//!   CCfits::VTshort, CCfits::VTuint, CCfits::VTint, CCfits::VTulong,
174//!   CCfits::VTlong, CCfits::VTlonglong, CCfits::VTfloat,
175//!   CCfits::VTdouble, CCfits::VTcomplex, CCfits::VTdblcomplex
176//!
177//! @param offsets
178//!         a vector containing the offsets to the columns (in bytes)
179//! @param targetFile
180//!         the ofstream where to write to
181//! @param fitsBuffer
182//!         the memory were the row has been loaded by cfitsio
183//
184int FitsDumper::WriteRow(ostream &out, const vector<CCfits::Column*> &list, const vector<int> &offsets, unsigned char* fitsBuffer) const
185{
186    int cnt = 0;
187
188    for (vector<CCfits::Column*>::const_iterator it=list.begin(); it!=list.end(); it++)
189    {
190        const CCfits::Column *col = *it;
191
192        // CCfits starts counting at 1 not 0
193        const int offset = offsets[col->index()-1];
194
195        // Get the pointer to the array which we are supposed to print
196        const unsigned char *ptr = fitsBuffer + offset;
197
198        // Loop over all array entries
199        for (int width=0; width<col->width(); width++)
200        {
201            switch (col->type())
202            {
203            case CCfits::Tbyte:     out << PtrToValue<uint8_t> (ptr); break;
204            case CCfits::Tushort:   out << PtrToValue<uint16_t>(ptr); break;
205            case CCfits::Tuint:
206            case CCfits::Tulong:    out << PtrToValue<uint32_t>(ptr); break;
207            case CCfits::Tshort:    out << PtrToValue<int16_t> (ptr); break;
208            case CCfits::Tint:
209            case CCfits::Tlong:     out << PtrToValue<int32_t> (ptr); break;
210            case CCfits::Tlonglong: out << PtrToValue<int64_t> (ptr); break;
211            case CCfits::Tfloat:    out << PtrToValue<float>   (ptr); break;
212            case CCfits::Tdouble:   out << PtrToValue<double>  (ptr); break;
213
214            default:
215                cerr << "Data type not implemented yet." << endl;
216                return 0;
217            }
218
219            out << " ";
220            cnt++;
221        }
222    }
223
224    if (cnt>0)
225        out << endl;
226
227    if (out.fail())
228    {
229        cerr << "ERROR - writing output: " << strerror(errno) << endl;
230        return -1;
231    }
232
233    return cnt;
234}
235
236// --------------------------------------------------------------------------
237//
238//! Calculates the required buffer size for reading one row of the current
239//! table. Also calculates the offsets to all the columns
240//
241vector<int> FitsDumper::CalculateOffsets() const
242{
243    map<int,int> sizes;
244
245    for (map<string, CCfits::Column*>::const_iterator it=fColMap.begin();
246         it!=fColMap.end(); it++)
247    {
248        const int &width = it->second->width();
249        const int &idx   = it->second->index();
250
251        const int size = ValueTypeToSize(it->second->type());
252        if (size==0)
253        {
254            cerr << "Data type " << (int)it->second->type() << " not implemented yet." << endl;
255            return vector<int>();
256        }
257
258        sizes[idx] = size*width;
259    }
260
261    //calculate the offsets in the vector.
262    vector<int> result(1, 0);
263
264    int size = 0;
265    int idx  = 0;
266
267    for (map<int,int>::const_iterator it=sizes.begin(); it!=sizes.end(); it++)
268    {
269        size += it->second;
270        result.push_back(size);
271
272        if (it->first == ++idx)
273            continue;
274
275        cerr << "Expected index " << idx << ", but found " << it->first << endl;
276        return vector<int>();
277    }
278
279    return result;
280}
281
282// --------------------------------------------------------------------------
283//
284//! Loads the fits file based on the current parameters
285//
286bool FitsDumper::OpenFile(const string &filename)
287{
288    if (fFile)
289        delete fFile;
290
291    ostringstream str;
292    try
293    {
294        fFile = new CCfits::FITS(filename);
295    }
296    catch (CCfits::FitsException e)
297    {
298        cerr << "Could not open FITS file " << filename << " reason: " << e.message() << endl;
299        return false;
300    }
301
302    return true;
303}
304
305bool FitsDumper::OpenTable(const string &tablename)
306{
307    if (!fFile)
308    {
309        cerr << "No file open." << endl;
310        return false;
311    }
312
313    const multimap< string, CCfits::ExtHDU * > extMap = fFile->extension();
314    if (extMap.find(tablename) == extMap.end())
315    {
316        cerr << "Table '" << tablename << "' not found." << endl;
317        return false;
318    }
319
320    fTable  = dynamic_cast<CCfits::Table*>(extMap.find(tablename)->second);
321    if (!fTable)
322    {
323        cerr << "Object '" << tablename << "' returned not a CCfits::Table." << endl;
324        return false;
325    }
326
327    fColMap = fTable->column();
328
329    fTable->makeThisCurrent();
330
331    fTable->getComments();
332    fTable->getHistory();
333    fTable->readAllKeys();
334
335    return true;
336}
337
338class MyColumn : public CCfits::Column
339{
340public:
341    const string &comment() const { return CCfits::Column::comment(); }
342};
343
344void FitsDumper::List()
345{
346    if (!fFile)
347    {
348        cerr << "No file open." << endl;
349        return;
350    }
351
352    cout << "\nFile: " << fFile->name() << "\n";
353
354    const multimap< string, CCfits::ExtHDU * > extMap = fFile->extension();
355    for (std::multimap<string, CCfits::ExtHDU*>::const_iterator it=extMap.begin(); it != extMap.end(); it++)
356    {
357
358        CCfits::Table *table = dynamic_cast<CCfits::Table*>(extMap.find(it->first)->second);
359
360        cout << " " << it->first << " [" << table->rows() << "]\n";
361
362        const map<string, CCfits::Column*> &cols = table->column();
363
364        for (map<string, CCfits::Column*>::const_iterator ic=cols.begin();
365             ic != cols.end(); ic++)
366        {
367            const MyColumn *col = static_cast<MyColumn*>(ic->second);
368
369            cout << "   " << col->name() << "[" << col->width() << "] * " << col->scale() << " (" << col->unit() << ":" <<  ValueTypeToStr(col->type())<<  ") " << col->comment() << "\n";
370            /*
371             inline size_t Column::repeat () const
372             inline bool   Column::varLength () const
373             inline double Column::zero () const
374             inline const  String& Column::display () const
375             inline const  String& Column::dimen () const
376             inline const  String& Column::TBCOL ()
377             inline const  String& Column::TTYPE ()
378             inline const  String& Column::TFORM ()
379             inline const  String& Column::TDISP ()
380             inline const  String& Column::TZERO ()
381             inline const  String& Column::TDIM  ()
382             inline const  String& Column::TNULL ()
383             inline const  String& Column::TLMIN ()
384             inline const  String& Column::TLMAX ()
385             inline const  String& Column::TDMAX ()
386             inline const  String& Column::TDMIN ()
387            */
388        }
389        cout << '\n';
390    }
391    cout << flush;
392}
393
394class MyKeyword : public CCfits::Keyword
395{
396public:
397    CCfits::ValueType keytype() const { return CCfits::Keyword::keytype(); }
398};
399
400void FitsDumper::ListKeywords(ostream &out)
401{
402    map<string,CCfits::Keyword*> keys = fTable->keyWord();
403    for (map<string,CCfits::Keyword*>::const_iterator it=keys.begin();
404         it!=keys.end(); it++)
405    {
406        string str;
407        double d;
408
409        const MyKeyword *kw = static_cast<MyKeyword*>(it->second);
410        kw->keytype();
411        out << "## " << setw(8) << kw->name() << "='";
412        if (kw->keytype()==16)
413            out << kw->value(str);
414        if (kw->keytype()==82)
415            out << kw->value(d);
416        out << "' \t(" << kw->comment() << ")" << endl;
417    }
418}
419
420void FitsDumper::ListHeader()
421{
422    if (!fTable)
423    {
424        cerr << "No table open." << endl;
425        return;
426    }
427
428    cout << "\nTable: " << fTable->name() << " (rows=" << fTable->rows() << ")\n";
429    if (!fTable->comment().empty())
430        cout << "Comment: \t" << fTable->comment() << '\n';
431    if (!fTable->history().empty())
432        cout << "History: \t" << fTable->history() << '\n';
433
434    ListKeywords(cout);
435    cout << endl;
436}
437
438
439// --------------------------------------------------------------------------
440//
441//! Perform the actual dump, based on the current parameters
442//
443bool FitsDumper::Dump(const string &filename, const vector<string> &list, int precision)
444{
445    for (vector<string>::const_iterator it=list.begin(); it!=list.end(); it++)
446        if (fColMap.find(*it) == fColMap.end())
447        {
448            cerr << "WARNING - Column '" << *it << "' not found in table." << endl;
449            return false;
450        }
451
452    // FIXME: Maybe do this when opening a table?
453    const vector<int> offsets = CalculateOffsets();
454    if (offsets.size()==0)
455        return false;
456
457    ofstream out(filename=="-"?"/dev/stdout":filename);
458    if (!out)
459    {
460        cerr << "Cannot open file " << filename << ": " << strerror(errno) << endl;
461        return false;
462    }
463
464    out.precision(precision);
465
466    out << "## --------------------------------------------------------------------------\n";
467    if (filename!="-")
468        out << "## File:    \t" << filename << '\n';
469    out << "## Table:   \t" << fTable->name() << '\n';
470    if (!fTable->comment().empty())
471        out << "## Comment: \t" << fTable->comment() << '\n';
472    if (!fTable->history().empty())
473        out << "## History: \t" << fTable->history() << '\n';
474    out << "## NumRows: \t" << fTable->rows() << '\n';
475    out << "## --------------------------------------------------------------------------\n";
476    ListKeywords(out);
477    out << "## --------------------------------------------------------------------------\n";
478    out << "#\n";
479    for (vector<string>::const_iterator it=list.begin(); it!=list.end(); it++)
480    {
481        const MyColumn *col = static_cast<MyColumn*>(fTable->column()[*it]);
482        out << "# " << col->name() << "[" << col->width() << "]: " << col->unit();
483        if (!col->comment().empty())
484            out << " (" <<col->comment() << ")";
485        out << '\n';
486    }
487    out << "#" << endl;
488
489
490    // Loop over all columns in our list of requested columns
491    vector<CCfits::Column*> columns;
492    for (vector<string>::const_iterator it=list.begin(); it!=list.end(); it++)
493        columns.push_back(fColMap.find(*it)->second);
494
495
496    const int size = offsets[offsets.size()-1];
497    unsigned char* fitsBuffer = new unsigned char[size];
498
499    int status = 0;
500    for (int i=1; i<=fTable->rows(); i++)
501    {
502        fits_read_tblbytes(fFile->fitsPointer(), i, 1, size, fitsBuffer, &status);
503        if (status)
504        {
505            cerr << "An error occurred while reading fits row #" << i << " error code: " << status << endl;
506            break;
507        }
508        if (WriteRow(out, columns, offsets, fitsBuffer)<0)
509        {
510            status=1;
511            break;
512        }
513    }
514    delete[] fitsBuffer;
515
516    return status==0;
517}
518
519// --------------------------------------------------------------------------
520//
521//! Perform the actual dump, based on the current parameters
522//
523/*
524bool FitsDumper::Plot(const vector<string> &list)
525{
526   for (vector<string>::const_iterator it=list.begin(); it!=list.end(); it++)
527        if (fColMap.find(*it) == fColMap.end())
528        {
529            cerr << "WARNING - Column '" << *it << "' not found in table." << endl;
530            return false;
531        }
532
533    // FIXME: Maybe do this when opening a table?
534    const vector<int> offsets = CalculateOffsets();
535    if (offsets.size()==0)
536        return false;
537
538    // Loop over all columns in our list of requested columns
539    const CCfits::Column *col[3] =
540    {
541        list.size()>0 ? fColMap.find(list[0])->second : 0,
542        list.size()>1 ? fColMap.find(list[1])->second : 0,
543        list.size()>2 ? fColMap.find(list[2])->second : 0
544    };
545
546    const int size = offsets[offsets.size()-1];
547    unsigned char* fitsBuffer = new unsigned char[size];
548
549    const int idx = 0;
550
551    // CCfits starts counting at 1 not 0
552    const size_t pos[3] =
553    {
554        col[0] ? offsets[col[0]->index()-1] + idx*col[0]->width()*ValueTypeToSize(col[0]->type()) : 0,
555        col[1] ? offsets[col[1]->index()-1] + idx*col[1]->width()*ValueTypeToSize(col[1]->type()) : 0,
556        col[2] ? offsets[col[2]->index()-1] + idx*col[2]->width()*ValueTypeToSize(col[2]->type()) : 0
557    };
558
559    int status = 0;
560    for (int i=1; i<=fTable->rows(); i++)
561    {
562        fits_read_tblbytes(fFile->fitsPointer(), i, 1, size, fitsBuffer, &status);
563        if (status)
564        {
565            cerr << "An error occurred while reading fits row #" << i << " error code: " << status << endl;
566            break;
567        }
568
569        try
570        {
571            const double x[3] =
572            {
573                col[0] ? PtrToDouble(fitsBuffer+pos[0], col[0]->type()) : 0,
574                col[1] ? PtrToDouble(fitsBuffer+pos[1], col[1]->type()) : 0,
575                col[2] ? PtrToDouble(fitsBuffer+pos[2], col[2]->type()) : 0
576            };
577        }
578        catch (const runtime_error &e)
579        {
580            cerr << e.what() << endl;
581            status=1;
582            break;
583        }
584    }
585    delete[] fitsBuffer;
586
587    return status==0;
588}
589*/
590
591// --------------------------------------------------------------------------
592//
593//! Retrieves the configuration parameters
594//! @param conf
595//!             the configuration object
596//
597int FitsDumper::ExecConfig(Configuration& conf)
598{
599    if (conf.Has("fitsfile"))
600    {
601        if (!OpenFile(conf.Get<string>("fitsfile")))
602            return -1;
603    }
604
605    if (conf.Get<bool>("list"))
606        List();
607
608    if (conf.Has("tablename"))
609    {
610        if (!OpenTable(conf.Get<string>("tablename")))
611            return -1;
612    }
613
614    if (conf.Get<bool>("header"))
615        ListHeader();
616
617    if (conf.Get<bool>("header") || conf.Get<bool>("list"))
618        return 1;
619
620    if (conf.Has("outfile"))
621    {
622        if (!Dump(conf.Get<string>("outfile"),
623                  conf.Get<vector<string>>("col"),
624                  conf.Get<int>("precision")))
625            return -1;
626    }
627
628    return 0;
629}
630
631void PrintUsage()
632{
633    cout <<
634        "fitsdump is a tool to dump data from a FITS table as ascii.\n"
635        "\n"
636        "Usage: fitsdump [OPTIONS] fitsfile col col ... \n"
637        "  or:  fitsdump [OPTIONS]\n";
638    cout << endl;
639}
640
641void PrintHelp()
642{
643    //
644}
645
646void SetupConfiguration(Configuration& conf)
647{
648    po::options_description configs("Fitsdump options");
649    configs.add_options()
650        ("fitsfile,f",  var<string>()
651#if BOOST_VERSION >= 104200
652         ->required()
653#endif
654                                                  , "Name of FITS file")
655        ("tablename,t", var<string>("DATA")
656#if BOOST_VERSION >= 104200
657         ->required()
658#endif
659                                                  , "Name of input table")
660        ("col,c",       vars<string>(),             "List of columns to dump")
661        ("outfile,o",   var<string>("/dev/stdout"), "Name of output file (-:/dev/stdout)")
662        ("precision,p", var<int>(20),               "Precision of ofstream")
663        ("list,l",      po_switch(),                "List all tables and columns in file")
664        ("header,h",    po_switch(),                "Dump header of given table")
665        ;
666
667    po::positional_options_description p;
668    p.add("fitsfile",   1); // The first positional options
669    p.add("col",       -1); // All others
670
671    conf.AddOptions(configs);
672    conf.SetArgumentPositions(p);
673}
674
675int main(int argc, const char** argv)
676{
677    Configuration conf(argv[0]);
678    conf.SetPrintUsage(PrintUsage);
679    SetupConfiguration(conf);
680
681    po::variables_map vm;
682    try
683    {
684        vm = conf.Parse(argc, argv);
685    }
686#if BOOST_VERSION > 104000
687    catch (po::multiple_occurrences &e)
688    {
689        cerr << "Program options invalid due to: " << e.what() << " of option '" << e.get_option_name() << "'." << endl;
690        return -1;
691    }
692#endif
693    catch (exception& e)
694    {
695        cerr << "Program options invalid due to: " << e.what() << endl;
696        return -1;
697    }
698
699    if (conf.HasVersion() || conf.HasPrint())
700        return -1;
701
702    if (conf.HasHelp())
703    {
704        PrintHelp();
705        return -1;
706    }
707
708    FitsDumper loader;
709    return loader.ExecConfig(conf);
710}
Note: See TracBrowser for help on using the repository browser.