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

Last change on this file since 12064 was 12060, checked in by tbretz, 13 years ago
Fixed includes.
File size: 34.1 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
15//#define PLOTTING_PLEASE
16
17#ifdef PLOTTING_PLEASE
18#include <QPen>
19#include <QtGui>
20#include <QApplication>
21
22#include <qwt_plot.h>
23#include <qwt_plot_grid.h>
24#include <qwt_plot_curve.h>
25#include <qwt_plot_zoomer.h>
26#include <qwt_legend.h>
27#include <qwt_scale_draw.h>
28
29#include "Time.h"
30#endif
31
32using namespace std;
33
34class MyColumn : public CCfits::Column
35{
36public:
37 const string &comment() const { return CCfits::Column::comment(); }
38 long width() const
39 {//the width() method seems to be buggy (or empty ?) as it always return 1... redo it ourselves.
40 string inter = format();
41 inter = inter.substr(0, inter.size()-1);
42 return atoi(inter.c_str());
43 }
44};
45
46#ifdef PLOTTING_PLEASE
47class TimeScaleDraw: public QwtScaleDraw
48{
49public:
50 virtual QwtText label(double v) const
51 {
52 Time t(v);
53 string time = t.GetAsStr("%H:%M:%S%F");
54 while (time[time.size()-1] == '0' && time.size() > 2)
55 {
56 time = time.substr(0, time.size()-1);
57 }
58 return QwtText(time.c_str());
59 }
60};
61#endif
62
63class FitsDumper
64{
65public:
66 FitsDumper();
67 ~FitsDumper();
68
69private:
70
71 CCfits::FITS *fFile; /// FITS pointer
72 CCfits::Table *fTable; /// Table pointer
73 map<string, MyColumn*> fColMap; /// map between the column names and their CCfits objects
74
75 // Convert CCfits::ValueType into a human readable string
76 string ValueTypeToStr(CCfits::ValueType type) const;
77
78 // Convert CCfits::ValueType into a number of associated bytes
79 int ValueTypeToSize(CCfits::ValueType type) const;
80
81 /// Calculate the buffer size required to read a row of the fits table, as well as the offsets to each column
82 vector<int> CalculateOffsets() const;
83
84 template<class T>
85 T PtrToValue(const unsigned char* &ptr) const;
86// template<class T>
87// double PtrToDouble(const unsigned char *ptr) const;
88// double PtrToDouble(const unsigned char *ptr, CCfits::ValueType type) const;
89
90 /// Write a single row of the selected data
91 int WriteRow(ostream &, const vector<MyColumn*> &, const vector<int> &, unsigned char *, const vector<pair<int, int> >&) const;
92
93 bool OpenFile(const string &); /// Open a file
94 bool OpenTable(const string &); /// Open a table
95
96 /// Lists all columns of an open file
97 void List();
98 void ListHeader();
99 void ListKeywords(ostream &);
100
101 bool separateColumnsFromRanges(const vector<string>& list,
102 vector<pair<int, int> >& ranges,
103 vector<string>& listNamesOnly);
104 /// Perform the dumping, based on the current dump list
105 bool Dump(const string &, const vector<string> &list, int);
106 ///Display the selected columns values VS time
107 int doCurvesDisplay( const vector<string> &list, const string& tableName);
108// bool Plot(const vector<string> &list);
109
110public:
111 ///Configures the fitsLoader from the config file and/or command arguments.
112 int ExecConfig(Configuration& conf);
113};
114
115// --------------------------------------------------------------------------
116//
117//! Constructor
118//! @param out
119//! the ostream where to redirect the outputs
120//
121FitsDumper::FitsDumper() : fFile(0), fTable(0)
122{
123}
124
125// --------------------------------------------------------------------------
126//
127//! Destructor
128//
129FitsDumper::~FitsDumper()
130{
131 if (fFile)
132 delete fFile;
133}
134
135
136string FitsDumper::ValueTypeToStr(CCfits::ValueType type) const
137{
138 switch (type)
139 {
140 case CCfits::Tbyte: return "uint8_t";
141 case CCfits::Tushort: return "uint16_t";
142 case CCfits::Tshort: return "int16_t";
143 case CCfits::Tuint: return "uint32_t";
144 case CCfits::Tint: return "int32_t";
145 case CCfits::Tulong: return "uint32_t";
146 case CCfits::Tlong: return "int32_t";
147 case CCfits::Tlonglong: return "int64_t";
148 case CCfits::Tfloat: return "float";
149 case CCfits::Tdouble: return "double";
150
151 default:
152 return "unknwown";
153 }
154}
155
156int FitsDumper::ValueTypeToSize(CCfits::ValueType type) const
157{
158 switch (type)
159 {
160 case CCfits::Tbyte: return sizeof(uint8_t);
161 case CCfits::Tushort:
162 case CCfits::Tshort: return sizeof(uint16_t);
163 case CCfits::Tuint:
164 case CCfits::Tint:
165 case CCfits::Tulong:
166 case CCfits::Tlong: return sizeof(uint32_t);
167 case CCfits::Tlonglong: return sizeof(uint64_t);
168 case CCfits::Tfloat: return sizeof(float);
169 case CCfits::Tdouble: return sizeof(double);
170
171 default:
172 return 0;
173 }
174}
175
176template<class T>
177T FitsDumper::PtrToValue(const unsigned char* &ptr) const
178{
179 T t;
180 reverse_copy(ptr, ptr+sizeof(T), reinterpret_cast<unsigned char*>(&t));
181 ptr += sizeof(T);
182
183 return t;
184}
185/*
186template<class T>
187double FitsDumper::PtrToDouble(const unsigned char *ptr) const
188{
189 T t;
190 reverse_copy(ptr, ptr+sizeof(T), reinterpret_cast<unsigned char*>(&t));
191 return t;
192}
193
194double FitsDumper::PtrToDouble(const unsigned char *ptr, CCfits::ValueType type) const
195{
196 switch (type)
197 {
198 case CCfits::Tbyte: return PtrToDouble<uint8_t> (ptr);
199 case CCfits::Tushort: return PtrToDouble<uint16_t>(ptr);
200 case CCfits::Tuint:
201 case CCfits::Tulong: return PtrToDouble<uint32_t>(ptr);
202 case CCfits::Tshort: return PtrToDouble<int16_t> (ptr);
203 case CCfits::Tint:
204 case CCfits::Tlong: return PtrToDouble<int32_t> (ptr);
205 case CCfits::Tlonglong: return PtrToDouble<int64_t> (ptr);
206 case CCfits::Tfloat: return PtrToDouble<float> (ptr);
207 case CCfits::Tdouble: return PtrToDouble<double> (ptr);
208
209 default:
210 throw runtime_error("Data type not implemented yet.");
211 }
212}
213*/
214
215// --------------------------------------------------------------------------
216//
217//! Writes a single row of the selected FITS data to the output file.
218//!
219//! Data type \b not yet implemented:
220//! CCfits::Tnull, CCfits::Tbit, CCfits::Tlogical, CCfits::Tstring,
221//! CCfits::Tcomplex, CCfits::Tdblcomplex, CCfits::VTbit,
222//! CCfits::VTbyte, CCfits::VTlogical, CCfits::VTushort,
223//! CCfits::VTshort, CCfits::VTuint, CCfits::VTint, CCfits::VTulong,
224//! CCfits::VTlong, CCfits::VTlonglong, CCfits::VTfloat,
225//! CCfits::VTdouble, CCfits::VTcomplex, CCfits::VTdblcomplex
226//!
227//! @param offsets
228//! a vector containing the offsets to the columns (in bytes)
229//! @param targetFile
230//! the ofstream where to write to
231//! @param fitsBuffer
232//! the memory were the row has been loaded by cfitsio
233//
234int FitsDumper::WriteRow(ostream &out, const vector<MyColumn*> &list, const vector<int> &offsets, unsigned char* fitsBuffer, const vector<pair<int, int> >& ranges) const
235{
236 int cnt = 0;
237 vector<pair<int, int> >::const_iterator jt = ranges.begin();
238 for (vector<MyColumn*>::const_iterator it=list.begin(); it!=list.end(); it++, jt++)
239 {
240 if (jt == ranges.end())
241 {
242 cout << "ERROR: END OF RANGE POINTER REACHED" << endl;
243 return false;
244 }
245 const MyColumn *col = *it;
246
247 // CCfits starts counting at 1 not 0
248 const int offset = offsets[col->index()-1];
249
250 // Get the pointer to the array which we are supposed to print
251 const unsigned char *ptr = fitsBuffer + offset;
252
253 // Loop over all array entries
254 int sizeToSkip = 0;
255 switch (col->type())
256 {
257 case CCfits::Tbyte: sizeToSkip = sizeof(uint8_t); break;
258 case CCfits::Tushort: sizeToSkip = sizeof(uint16_t); break;
259 case CCfits::Tuint:
260 case CCfits::Tulong: sizeToSkip = sizeof(uint32_t); break;
261 case CCfits::Tshort: sizeToSkip = sizeof(int16_t); break;
262 case CCfits::Tint:
263 case CCfits::Tlong: sizeToSkip = sizeof(int32_t); break;
264 case CCfits::Tlonglong: sizeToSkip = sizeof(int64_t); break;
265 case CCfits::Tfloat: sizeToSkip = sizeof(float); break;
266 case CCfits::Tdouble: sizeToSkip = sizeof(double); break;
267 default:
268 cerr << "Data type not implemented yet." << endl;
269 return 0;
270 }
271 ptr += sizeToSkip*jt->first;
272 for (int width=jt->first; width<jt->second; width++)
273 {
274 switch (col->type())
275 {
276 case CCfits::Tbyte: out << PtrToValue<uint8_t> (ptr); break;
277 case CCfits::Tushort: out << PtrToValue<uint16_t>(ptr); break;
278 case CCfits::Tuint:
279 case CCfits::Tulong: out << PtrToValue<uint32_t>(ptr); break;
280 case CCfits::Tshort: out << PtrToValue<int16_t> (ptr); break;
281 case CCfits::Tint:
282 case CCfits::Tlong: out << PtrToValue<int32_t> (ptr); break;
283 case CCfits::Tlonglong: out << PtrToValue<int64_t> (ptr); break;
284 case CCfits::Tfloat: out << PtrToValue<float> (ptr); break;
285 case CCfits::Tdouble: out << PtrToValue<double> (ptr); break;
286
287 default:
288 cerr << "Data type not implemented yet." << endl;
289 return 0;
290 }
291
292 out << " ";
293 cnt++;
294 }
295 }
296
297 if (cnt>0)
298 out << endl;
299
300 if (out.fail())
301 {
302 cerr << "ERROR - writing output: " << strerror(errno) << endl;
303 return -1;
304 }
305
306 return cnt;
307}
308
309// --------------------------------------------------------------------------
310//
311//! Calculates the required buffer size for reading one row of the current
312//! table. Also calculates the offsets to all the columns
313//
314vector<int> FitsDumper::CalculateOffsets() const
315{
316 map<int,int> sizes;
317
318 for (map<string, MyColumn*>::const_iterator it=fColMap.begin();
319 it!=fColMap.end(); it++)
320 {
321 const int &width = it->second->width();
322 const int &idx = it->second->index();
323
324 const int size = ValueTypeToSize(it->second->type());
325 if (size==0)
326 {
327 cerr << "Data type " << (int)it->second->type() << " not implemented yet." << endl;
328 return vector<int>();
329 }
330
331 sizes[idx] = size*width;
332 }
333
334 //calculate the offsets in the vector.
335 vector<int> result(1, 0);
336
337 int size = 0;
338 int idx = 0;
339
340 for (map<int,int>::const_iterator it=sizes.begin(); it!=sizes.end(); it++)
341 {
342 size += it->second;
343 result.push_back(size);
344
345 if (it->first == ++idx)
346 continue;
347
348 cerr << "Expected index " << idx << ", but found " << it->first << endl;
349 return vector<int>();
350 }
351
352 return result;
353}
354
355// --------------------------------------------------------------------------
356//
357//! Loads the fits file based on the current parameters
358//
359bool FitsDumper::OpenFile(const string &filename)
360{
361 if (fFile)
362 delete fFile;
363
364 ostringstream str;
365 try
366 {
367 fFile = new CCfits::FITS(filename);
368 }
369 catch (CCfits::FitsException e)
370 {
371 cerr << "Could not open FITS file " << filename << " reason: " << e.message() << endl;
372 return false;
373 }
374
375 return true;
376}
377
378bool FitsDumper::OpenTable(const string &tablename)
379{
380 if (!fFile)
381 {
382 cerr << "No file open." << endl;
383 return false;
384 }
385
386 const multimap< string, CCfits::ExtHDU * > extMap = fFile->extension();
387 if (extMap.find(tablename) == extMap.end())
388 {
389 cerr << "Table '" << tablename << "' not found." << endl;
390 return false;
391 }
392
393 fTable = dynamic_cast<CCfits::Table*>(extMap.find(tablename)->second);
394 if (!fTable)
395 {
396 cerr << "Object '" << tablename << "' returned not a CCfits::Table." << endl;
397 return false;
398 }
399
400 map<string, CCfits::Column*>& tCols = fTable->column();
401 for (map<string, CCfits::Column*>::const_iterator it = tCols.begin(); it != tCols.end(); it++)
402 {
403 fColMap.insert(make_pair(it->first, static_cast<MyColumn*>(it->second)));
404 }
405// fColMap = fTable->column();
406
407 fTable->makeThisCurrent();
408
409 fTable->getComments();
410 fTable->getHistory();
411 fTable->readAllKeys();
412
413 return true;
414}
415
416
417void FitsDumper::List()
418{
419 if (!fFile)
420 {
421 cerr << "No file open." << endl;
422 return;
423 }
424
425 cout << "\nFile: " << fFile->name() << "\n";
426
427 const multimap< string, CCfits::ExtHDU * > extMap = fFile->extension();
428 for (std::multimap<string, CCfits::ExtHDU*>::const_iterator it=extMap.begin(); it != extMap.end(); it++)
429 {
430
431 CCfits::Table *table = dynamic_cast<CCfits::Table*>(extMap.find(it->first)->second);
432
433 table->makeThisCurrent();
434 table->readData();
435
436 cout << " " << it->first << " [" << table->rows() << "]\n";
437
438 const map<string, CCfits::Column*> &cols = table->column();
439
440// for (map<string, CCfits::Column*>::const_iterator id = cols.begin(); ib != cols.end(); ib++)
441// {
442// TFORM
443// }
444
445 for (map<string, CCfits::Column*>::const_iterator ic=cols.begin();
446 ic != cols.end(); ic++)
447 {
448 const MyColumn *col = static_cast<MyColumn*>(ic->second);
449
450 cout << " " << col->name() << "[" << col->width() << "] * " << col->scale() << " (" << col->unit() << ":" << ValueTypeToStr(col->type())<< ") " << col->comment() << "\n";
451 /*
452 inline size_t Column::repeat () const
453 inline bool Column::varLength () const
454 inline double Column::zero () const
455 inline const String& Column::display () const
456 inline const String& Column::dimen () const
457 inline const String& Column::TBCOL ()
458 inline const String& Column::TTYPE ()
459 inline const String& Column::TFORM ()
460 inline const String& Column::TDISP ()
461 inline const String& Column::TZERO ()
462 inline const String& Column::TDIM ()
463 inline const String& Column::TNULL ()
464 inline const String& Column::TLMIN ()
465 inline const String& Column::TLMAX ()
466 inline const String& Column::TDMAX ()
467 inline const String& Column::TDMIN ()
468 */
469 }
470 cout << '\n';
471 }
472 cout << flush;
473}
474
475class MyKeyword : public CCfits::Keyword
476{
477public:
478 CCfits::ValueType keytype() const { return CCfits::Keyword::keytype(); }
479};
480
481void FitsDumper::ListKeywords(ostream &out)
482{
483 map<string,CCfits::Keyword*> keys = fTable->keyWord();
484 for (map<string,CCfits::Keyword*>::const_iterator it=keys.begin();
485 it!=keys.end(); it++)
486 {
487 string str;
488 double d;
489 int l;
490
491 const MyKeyword *kw = static_cast<MyKeyword*>(it->second);
492 kw->keytype();
493 out << "## " << setw(8) << kw->name() << " = " << setw(10);
494 if (kw->keytype()==16)
495 out << ("'"+kw->value(str)+"'");
496 if (kw->keytype()==31)
497 out << kw->value(l);
498 if (kw->keytype()==82)
499 out << kw->value(d);
500 out << " / " << kw->comment() << endl;
501 }
502}
503
504void FitsDumper::ListHeader()
505{
506 if (!fTable)
507 {
508 cerr << "No table open." << endl;
509 return;
510 }
511
512 cout << "\nTable: " << fTable->name() << " (rows=" << fTable->rows() << ")\n";
513 if (!fTable->comment().empty())
514 cout << "Comment: \t" << fTable->comment() << '\n';
515 if (!fTable->history().empty())
516 cout << "History: \t" << fTable->history() << '\n';
517
518 ListKeywords(cout);
519 cout << endl;
520}
521
522bool FitsDumper::separateColumnsFromRanges(const vector<string>& list,
523 vector<pair<int, int> >& ranges,
524 vector<string>& listNamesOnly)
525{
526 for (vector<string>::const_iterator it=list.begin(); it!=list.end(); it++)
527 {
528 string columnNameOnly = *it;
529 unsigned long bracketIndex0 = columnNameOnly.find_first_of('[');
530 unsigned long bracketIndex1 = columnNameOnly.find_first_of(']');
531 unsigned long colonIndex = columnNameOnly.find_first_of(':');
532// cout << bracketIndex0 << " " << bracketIndex1 << " " << colonIndex << endl;
533 int columnStart = -1;
534 int columnEnd = -1;
535 if (bracketIndex0 != string::npos)
536 {//there is a range given. Extract the range
537 if (colonIndex != string::npos)
538 {//we have a range here
539 columnStart = atoi(columnNameOnly.substr(bracketIndex0+1, colonIndex-(bracketIndex0+1)).c_str());
540 columnEnd = atoi(columnNameOnly.substr(colonIndex+1, bracketIndex1-(colonIndex+1)).c_str());
541 columnEnd++;
542 }
543 else
544 {//only a single index there
545 columnStart = atoi(columnNameOnly.substr(bracketIndex0+1, bracketIndex1 - (bracketIndex0+1)).c_str());
546 columnEnd = columnStart+1;
547// cout << "Cstart " << columnStart << " end: " << columnEnd << endl;
548 }
549 columnNameOnly = columnNameOnly.substr(0, bracketIndex0);
550 }
551
552 if (fColMap.find(columnNameOnly) == fColMap.end())
553 {
554 cerr << "ERROR - Column '" << columnNameOnly << "' not found in table." << endl;
555 return false;
556 }
557// cout << "The column name is: " << columnNameOnly << endl;
558 MyColumn *cCol = static_cast<MyColumn*>(fColMap.find(columnNameOnly)->second);
559 if (bracketIndex0 == string::npos)
560 {//no range given: use the full range
561 ranges.push_back(make_pair(0, cCol->width()));
562 columnStart = 0;
563 columnEnd = 1;
564 }
565 else
566 {//use the range extracted earlier
567 if (columnStart < 0)
568 {
569 cerr << "ERROR - Start range for column " << columnNameOnly << " is less than zero (" << columnStart << "). Aborting" << endl;
570 return false;
571 }
572 if (columnEnd > cCol->width())
573 {
574 cerr << "ERROR - End range for column " << columnNameOnly << " is greater than the last element (" << cCol->width() << " vs " << columnEnd << "). Aborting" << endl;
575 return false;
576 }
577 ranges.push_back(make_pair(columnStart, columnEnd));
578 }
579// cout << "Will be exporting from " << columnStart << " to " << columnEnd-1 << " for column " << columnNameOnly << endl;
580 listNamesOnly.push_back(columnNameOnly);
581 }
582 return true;
583}
584// --------------------------------------------------------------------------
585//
586//! Perform the actual dump, based on the current parameters
587//
588bool FitsDumper::Dump(const string &filename, const vector<string> &list, int precision)
589{
590
591 //first of all, let's separate the columns from their ranges and check that the requested columns are indeed part of the file
592 vector<pair<int, int> > ranges;
593 vector<string> listNamesOnly;
594
595 if (!separateColumnsFromRanges(list, ranges, listNamesOnly))
596 {
597 cerr << "Something went wrong while extracting the columns names from parameters. Aborting" << endl;
598 return false;
599 }
600
601 // FIXME: Maybe do this when opening a table?
602 const vector<int> offsets = CalculateOffsets();
603 if (offsets.size()==0)
604 return false;
605
606 ofstream out(filename=="-"?"/dev/stdout":filename);
607 if (!out)
608 {
609 cerr << "Cannot open file " << filename << ": " << strerror(errno) << endl;
610 return false;
611 }
612
613 out.precision(precision);
614
615 out << "## --------------------------------------------------------------------------\n";
616 if (filename!="-")
617 out << "## File: \t" << filename << '\n';
618 out << "## Table: \t" << fTable->name() << '\n';
619 if (!fTable->comment().empty())
620 out << "## Comment: \t" << fTable->comment() << '\n';
621 if (!fTable->history().empty())
622 out << "## History: \t" << fTable->history() << '\n';
623 out << "## NumRows: \t" << fTable->rows() << '\n';
624 out << "## --------------------------------------------------------------------------\n";
625 ListKeywords(out);
626 out << "## --------------------------------------------------------------------------\n";
627 out << "#\n";
628 //vector<pair<int, int> >::const_iterator rangesIt = ranges.begin();
629 //the above is soooo yesterday ;) let's try the new auto feature of c++0x
630 auto rangesIt = ranges.begin();
631 for (vector<string>::const_iterator it=listNamesOnly.begin(); it!=listNamesOnly.end(); it++, rangesIt++)
632 {
633 const MyColumn *col = static_cast<MyColumn*>(fTable->column()[*it]);
634 if (rangesIt->first != 0 || rangesIt->second != col->width())
635 {
636 out << "#";
637 for (int i=rangesIt->first; i<rangesIt->second; i++)
638 out << " " << col->name() << "[" << i << "]";
639 out << ": " << col->unit();
640 }
641 else
642 out << "# " << col->name() << "[" << col->width() << "]: " << col->unit();
643
644 if (!col->comment().empty())
645 out << " (" <<col->comment() << ")";
646 out << '\n';
647 }
648 out << "#" << endl;
649
650
651 // Loop over all columns in our list of requested columns
652 vector<MyColumn*> columns;
653 for (vector<string>::const_iterator it=listNamesOnly.begin(); it!=listNamesOnly.end(); it++)
654 {
655 MyColumn *cCol = static_cast<MyColumn*>(fColMap.find(*it)->second);
656 columns.push_back(cCol);
657 }
658 const int size = offsets[offsets.size()-1];
659 unsigned char* fitsBuffer = new unsigned char[size];
660
661 int status = 0;
662 for (int i=1; i<=fTable->rows(); i++)
663 {
664 fits_read_tblbytes(fFile->fitsPointer(), i, 1, size, fitsBuffer, &status);
665 if (status)
666 {
667 cerr << "An error occurred while reading fits row #" << i << " error code: " << status << endl;
668 break;
669 }
670 if (WriteRow(out, columns, offsets, fitsBuffer, ranges)<0)
671 {
672 status=1;
673 break;
674 }
675 }
676 delete[] fitsBuffer;
677
678 return status==0;
679}
680
681// --------------------------------------------------------------------------
682//
683//! Perform the actual dump, based on the current parameters
684//
685/*
686bool FitsDumper::Plot(const vector<string> &list)
687{
688 for (vector<string>::const_iterator it=list.begin(); it!=list.end(); it++)
689 if (fColMap.find(*it) == fColMap.end())
690 {
691 cerr << "WARNING - Column '" << *it << "' not found in table." << endl;
692 return false;
693 }
694
695 // FIXME: Maybe do this when opening a table?
696 const vector<int> offsets = CalculateOffsets();
697 if (offsets.size()==0)
698 return false;
699
700 // Loop over all columns in our list of requested columns
701 const CCfits::Column *col[3] =
702 {
703 list.size()>0 ? fColMap.find(list[0])->second : 0,
704 list.size()>1 ? fColMap.find(list[1])->second : 0,
705 list.size()>2 ? fColMap.find(list[2])->second : 0
706 };
707
708 const int size = offsets[offsets.size()-1];
709 unsigned char* fitsBuffer = new unsigned char[size];
710
711 const int idx = 0;
712
713 // CCfits starts counting at 1 not 0
714 const size_t pos[3] =
715 {
716 col[0] ? offsets[col[0]->index()-1] + idx*col[0]->width()*ValueTypeToSize(col[0]->type()) : 0,
717 col[1] ? offsets[col[1]->index()-1] + idx*col[1]->width()*ValueTypeToSize(col[1]->type()) : 0,
718 col[2] ? offsets[col[2]->index()-1] + idx*col[2]->width()*ValueTypeToSize(col[2]->type()) : 0
719 };
720
721 int status = 0;
722 for (int i=1; i<=fTable->rows(); i++)
723 {
724 fits_read_tblbytes(fFile->fitsPointer(), i, 1, size, fitsBuffer, &status);
725 if (status)
726 {
727 cerr << "An error occurred while reading fits row #" << i << " error code: " << status << endl;
728 break;
729 }
730
731 try
732 {
733 const double x[3] =
734 {
735 col[0] ? PtrToDouble(fitsBuffer+pos[0], col[0]->type()) : 0,
736 col[1] ? PtrToDouble(fitsBuffer+pos[1], col[1]->type()) : 0,
737 col[2] ? PtrToDouble(fitsBuffer+pos[2], col[2]->type()) : 0
738 };
739 }
740 catch (const runtime_error &e)
741 {
742 cerr << e.what() << endl;
743 status=1;
744 break;
745 }
746 }
747 delete[] fitsBuffer;
748
749 return status==0;
750}
751*/
752
753// --------------------------------------------------------------------------
754//
755//! Retrieves the configuration parameters
756//! @param conf
757//! the configuration object
758//
759int FitsDumper::ExecConfig(Configuration& conf)
760{
761 if (conf.Has("fitsfile"))
762 {
763 if (!OpenFile(conf.Get<string>("fitsfile")))
764 return -1;
765 }
766
767 if (conf.Get<bool>("list"))
768 List();
769
770 if (conf.Has("tablename"))
771 {
772 if (!OpenTable(conf.Get<string>("tablename")))
773 return -1;
774 }
775
776 if (conf.Get<bool>("graph"))
777 {
778 if (!conf.Has("col"))
779 {
780 cout << "Please specify the columns that should be dumped as arguments. Aborting" << endl;
781 return 0;
782 }
783 doCurvesDisplay(conf.Get<vector<string>>("col"),
784 conf.Get<string>("tablename"));
785 return 1;
786 }
787
788
789 if (conf.Get<bool>("header"))
790 ListHeader();
791
792 if (conf.Get<bool>("header") || conf.Get<bool>("list"))
793 return 1;
794
795 if (conf.Has("outfile"))
796 {
797 if (!conf.Has("col"))
798 {
799 cout << "Please specify the columns that should be dumped as arguments. Aborting" << endl;
800 return 0;
801 }
802 if (!Dump(conf.Get<string>("outfile"),
803 conf.Get<vector<string>>("col"),
804 conf.Get<int>("precision")))
805 return -1;
806 }
807
808
809 return 0;
810}
811
812void PrintUsage()
813{
814 cout <<
815 "fitsdump is a tool to dump data from a FITS table as ascii.\n"
816 "\n"
817 "Usage: fitsdump [OPTIONS] fitsfile col col ... \n"
818 " or: fitsdump [OPTIONS]\n";
819 cout << endl;
820}
821
822void PrintHelp()
823{
824 //
825}
826#ifdef PLOTTING_PLEASE
827int FitsDumper::doCurvesDisplay( const vector<string> &list, const string& tableName)
828{
829 //first of all, let's separate the columns from their ranges and check that the requested columns are indeed part of the file
830 vector<pair<int, int> > ranges;
831 vector<string> listNamesOnly;
832 if (!separateColumnsFromRanges(list, ranges, listNamesOnly))
833 {
834 cerr << "Something went wrong while extracting the columns names from parameters. Aborting" << endl;
835 return false;
836 }
837 vector<string> curvesNames;
838 stringstream str;
839 for (auto it=ranges.begin(), jt=listNamesOnly.begin(); it != ranges.end(); it++, jt++)
840 {
841 for (int i=it->first; i<it->second;i++)
842 {
843 str.str("");
844 str << *jt << "[" << i << "]";
845 curvesNames.push_back(str.str());
846 }
847 }
848 char* handle = new char[17];
849 sprintf(handle,"FitsDump Display");
850// Qt::HANDLE h = *handle;//NULL
851 int argc = 1;
852 char ** argv = &handle;
853 QApplication a(argc, argv);
854
855
856
857 QwtPlot* plot = new QwtPlot();
858 QwtPlotGrid* grid = new QwtPlotGrid;
859 grid->enableX(false);
860 grid->enableY(true);
861 grid->enableXMin(false);
862 grid->enableYMin(false);
863 grid->setMajPen(QPen(Qt::black, 0, Qt::DotLine));
864 grid->attach(plot);
865 plot->setAutoReplot(true);
866 string title = tableName;
867 plot->setAxisScaleDraw( QwtPlot::xBottom, new TimeScaleDraw());
868
869 QWidget window;
870 QHBoxLayout* layout = new QHBoxLayout(&window);
871 layout->setContentsMargins(0,0,0,0);
872 layout->addWidget(plot);
873
874 QwtPlotZoomer zoom(plot->canvas());
875 zoom.setRubberBandPen(QPen(Qt::gray, 2, Qt::DotLine));
876 zoom.setTrackerPen(QPen(Qt::gray));
877 int totalSize = 0;
878 for (unsigned int i=0;i<list.size();i++)
879 totalSize += ranges[i].second - ranges[i].first;
880 cout << "Total size: " << totalSize << endl;
881 vector<QwtPlotCurve*> curves(totalSize);
882 int ii=0;
883 for (auto it = curves.begin(), jt=curvesNames.begin(); it != curves.end(); it++, jt++)
884 {
885 *it = new QwtPlotCurve(jt->c_str());
886 switch (ii%6)
887 {
888 case 0:
889 (*it)->setPen(QColor(255,0,0));
890 break;
891 case 1:
892 (*it)->setPen(QColor(0,255,0));
893 break;
894 case 2:
895 (*it)->setPen(QColor(0,0,255));
896 break;
897 case 3:
898 (*it)->setPen(QColor(255,255,0));
899 break;
900 case 4:
901 (*it)->setPen(QColor(0,255,255));
902 break;
903 case 5:
904 (*it)->setPen(QColor(255,0,255));
905 break;
906 default:
907 (*it)->setPen(QColor(0,0,0));
908 };
909 ii++;
910 (*it)->setStyle(QwtPlotCurve::Lines);
911 (*it)->attach(plot);
912 }
913 plot->insertLegend(new QwtLegend(), QwtPlot::RightLegend);
914
915 const vector<int> offsets = CalculateOffsets();
916 if (offsets.size()==0)
917 return false;
918
919
920 // Loop over all columns in our list of requested columns
921 vector<MyColumn*> columns;
922 for (vector<string>::const_iterator it=listNamesOnly.begin(); it!=listNamesOnly.end(); it++)
923 {
924 MyColumn *cCol = static_cast<MyColumn*>(fColMap.find(*it)->second);
925 columns.push_back(cCol);
926 }
927 //add the time column to the given columns
928 MyColumn* timeCol = static_cast<MyColumn*>(fColMap.find("Time")->second);
929 if (!timeCol)
930 {
931 cerr << "Error: time column could not be found in given table. Aborting" << endl;
932 return false;
933 }
934 columns.push_back(timeCol);
935 ranges.push_back(make_pair(0,1));
936 /////
937 const int size = offsets[offsets.size()-1];
938 unsigned char* fitsBuffer = new unsigned char[size];
939
940// stringstream str;
941 str.str("");
942 int status = 0;
943
944 vector<double*> xValues(totalSize);
945 double* yValues;
946 cout.precision(10);
947 str.precision(20);
948 for (auto it=xValues.begin(); it!=xValues.end(); it++)
949 *it = new double[fTable->rows()];
950
951 yValues = new double[fTable->rows()];
952
953 for (int i=1; i<=fTable->rows(); i++)
954 {
955 fits_read_tblbytes(fFile->fitsPointer(), i, 1, size, fitsBuffer, &status);
956 if (status)
957 {
958 cerr << "An error occurred while reading fits row #" << i << " error code: " << status << endl;
959 break;
960 }
961 if (WriteRow(str, columns, offsets, fitsBuffer, ranges)<0)
962 {
963 status=1;
964 cerr << "An Error occured while reading the fits row " << i << endl;
965 return -1;
966 }
967// yValues[i-1] = i;
968 for (auto it=xValues.begin(); it!= xValues.end(); it++)
969 {
970 str >> (*it)[i-1];
971// cout << (*it)[i-1] << " ";
972 }
973 str >> yValues[i-1];
974 if (i==1)
975 {
976 Time t(yValues[0]);
977 title += " - " + t.GetAsStr("%Y-%m-%d");
978 plot->setTitle(title.c_str());
979 }
980// cout << yValues[i-1] << " ";
981// cout << endl;
982 }
983 //set the actual data.
984 auto jt = xValues.begin();
985 for (auto it=curves.begin(); it != curves.end(); it++, jt++)
986 (*it)->setRawData(yValues, *jt, fTable->rows());
987
988 QStack<QRectF> stack;
989 double minX, minY, maxX, maxY;
990 minX = minY = 1e10;
991 maxX = maxY = -1e10;
992 QRectF rect;
993 QPointF point;
994 for (auto it=curves.begin(); it!= curves.end(); it++)
995 {
996 rect = (*it)->boundingRect();
997 point = rect.bottomRight();
998 if (point.x() < minX) minX = point.x();
999 if (point.y() < minY) minY = point.y();
1000 if (point.x() > maxX) maxX = point.x();
1001 if (point.y() > maxY) maxY = point.y();
1002 point = rect.topLeft();
1003 if (point.x() < minX) minX = point.x();
1004 if (point.y() < minY) minY = point.y();
1005 if (point.x() > maxX) maxX = point.x();
1006 if (point.y() > maxY) maxY = point.y();
1007 }
1008 QPointF bottomRight(maxX, minY);
1009 QPointF topLeft(minX, maxY);
1010 QPointF center((bottomRight+topLeft)/2.f);
1011 stack.push(QRectF(topLeft + (topLeft-center)*(.5f),bottomRight + (bottomRight-center)*(.5f)));
1012 zoom.setZoomStack(stack);
1013
1014 delete[] fitsBuffer;
1015 window.resize(600, 400);
1016 window.show();
1017
1018 a.exec();
1019
1020 for (auto it = curves.begin(); it != curves.end(); it++)
1021 delete *it;
1022 for (auto it = xValues.begin(); it != xValues.end(); it++)
1023 delete[] *it;
1024 delete[] yValues;
1025 delete[] handle;
1026 return 0;
1027}
1028#else
1029int FitsDumper::doCurvesDisplay( const vector<string> &list, const string& tableName)
1030{
1031 cerr << "Sorry, but plotting features seem to have been disabled at compilation time." << endl;
1032 cerr << "Please recompile with PLOTTING_PLEASE defined and try again." << endl;
1033 return 0;
1034}
1035#endif
1036void SetupConfiguration(Configuration& conf)
1037{
1038 po::options_description configs("Fitsdump options");
1039 configs.add_options()
1040 ("fitsfile,f", var<string>()
1041#if BOOST_VERSION >= 104200
1042 ->required()
1043#endif
1044 , "Name of FITS file")
1045 ("tablename,t", var<string>("DATA")
1046#if BOOST_VERSION >= 104200
1047 ->required()
1048#endif
1049 , "Name of input table")
1050 ("col,c", vars<string>(), "List of columns to dump\narg is a list of columns, separated by a space.\nAdditionnally, a list of sub-columns can be added\ne.g. Data[3] will dump sub-column 3 of column Data\nData[3:4] will dump sub-columns 3 and 4\nOmitting this argument dump the entire column\nnota: all indices start at zero")
1051 ("outfile,o", var<string>("/dev/stdout"), "Name of output file (-:/dev/stdout)")
1052 ("precision,p", var<int>(20), "Precision of ofstream")
1053 ("list,l", po_switch(), "List all tables and columns in file")
1054 ("header,h", po_switch(), "Dump header of given table")
1055 ("graph,g", po_switch(), "Plot the columns instead of dumping them")
1056 ;
1057
1058 po::positional_options_description p;
1059 p.add("fitsfile", 1); // The first positional options
1060 p.add("col", -1); // All others
1061
1062 conf.AddOptions(configs);
1063 conf.SetArgumentPositions(p);
1064}
1065
1066int main(int argc, const char** argv)
1067{
1068 Configuration conf(argv[0]);
1069 conf.SetPrintUsage(PrintUsage);
1070 SetupConfiguration(conf);
1071
1072 if (!conf.DoParse(argc, argv, PrintHelp))
1073 return -1;
1074
1075 FitsDumper loader;
1076 return loader.ExecConfig(conf);
1077}
Note: See TracBrowser for help on using the repository browser.