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

Last change on this file since 12725 was 12725, checked in by lyard, 13 years ago
bugfix fitsdump
File size: 31.7 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 "externals/fits.h"
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
34
35#ifdef PLOTTING_PLEASE
36class TimeScaleDraw: public QwtScaleDraw
37{
38public:
39 virtual QwtText label(double v) const
40 {
41 Time t(v);
42 string time = t.GetAsStr("%H:%M:%S%F");
43 while (time[time.size()-1] == '0' && time.size() > 2)
44 {
45 time = time.substr(0, time.size()-1);
46 }
47 return QwtText(time.c_str());
48 }
49};
50#endif
51
52class FitsDumper
53{
54public:
55 FitsDumper();
56 ~FitsDumper();
57
58private:
59 fits* fFile;
60 bool fDotsPlease;
61 string fFilename;
62
63 fits::Table::Columns fColMap;
64 fits::Table::Keys fKeyMap;
65
66 // Convert CCfits::ValueType into a human readable string
67 string ValueTypeToStr(char type) const;
68
69 // Convert CCfits::ValueType into a number of associated bytes
70 int ValueTypeToSize(char type) const;
71
72 /// Calculate the buffer size required to read a row of the fits table, as well as the offsets to each column
73// vector<int> CalculateOffsets() const;
74
75 template<class T>
76 T PtrToValue(const unsigned char* &ptr) const;
77// template<class T>
78// double PtrToDouble(const unsigned char *ptr) const;
79// double PtrToDouble(const unsigned char *ptr, CCfits::ValueType type) const;
80
81 /// Write a single row of the selected data
82 // int WriteRow(ostream &, const vector<MyColumn*> &, const vector<int> &, unsigned char *, const vector<pair<int, int> >&) const;
83
84 bool OpenFile(const string &); /// Open a file
85 bool OpenTable(const string &); /// Open a table
86
87 /// Lists all columns of an open file
88 void List();
89 void ListHeader();
90 void ListKeywords(ostream &);
91
92 bool separateColumnsFromRanges(const vector<string>& list,
93 vector<pair<int, int> >& ranges,
94 vector<string>& listNamesOnly);
95 /// Perform the dumping, based on the current dump list
96 bool Dump(const string &, const vector<string> &list, int);
97 ///Display the selected columns values VS time
98#ifdef PLOTTING_PLEASE
99 int doCurvesDisplay( const vector<string> &list, const string& tableName);
100#endif
101 int doStatsPlease(const string &filename, const vector<string>& list, int precision);
102 // bool Plot(const vector<string> &list);
103
104public:
105 ///Configures the fitsLoader from the config file and/or command arguments.
106 int ExecConfig(Configuration& conf);
107};
108
109// --------------------------------------------------------------------------
110//
111//! Constructor
112//! @param out
113//! the ostream where to redirect the outputs
114//
115FitsDumper::FitsDumper() : fFile(0), fDotsPlease(false)
116{
117}
118
119// --------------------------------------------------------------------------
120//
121//! Destructor
122//
123FitsDumper::~FitsDumper()
124{
125 if (fFile)
126 delete fFile;
127}
128
129
130string FitsDumper::ValueTypeToStr(char type) const
131{
132 switch (type)
133 {
134 case 'L': return "bool(8)";
135 case 'B': return "byte(8)";
136 case 'I': return "short(16)";
137 case 'J': return "int(32)";
138 case 'K': return "int(64)";
139 case 'E': return "float(32)";
140 case 'D': return "double(64)";
141 default:
142 return "unknwown";
143 }
144}
145
146int FitsDumper::ValueTypeToSize(char type) const
147{
148 switch (type)
149 {
150 case 'L': return sizeof(uint8_t);
151 case 'B': return sizeof(int8_t);
152 case 'I': return sizeof(int16_t);
153 case 'J': return sizeof(int32_t);
154 case 'K': return sizeof(int64_t);
155 case 'E': return sizeof(float);
156 case 'D': return sizeof(double);
157 default:
158 return 0;
159 }
160}
161
162template<class T>
163T FitsDumper::PtrToValue(const unsigned char* &ptr) const
164{
165 T t;
166 reverse_copy(ptr, ptr+sizeof(T), reinterpret_cast<unsigned char*>(&t));
167 ptr += sizeof(T);
168
169 return t;
170}
171// --------------------------------------------------------------------------
172//
173//! Loads the fits file based on the current parameters
174//
175bool FitsDumper::OpenFile(const string &filename)
176{
177 if (fFile)
178 {
179 fFile->close();
180 delete fFile;
181 }
182
183 try {
184 fFile = new fits(filename);
185 }
186 catch (std::runtime_error e)
187 {
188 cout << "Something went wrong while trying to open " << filename;
189 cout << ": " << e.what() << " Aborting dump." << endl;
190 return false;
191 }
192 fFilename = filename;
193
194 const fits::Table::Columns& tCols = fFile->getColumns();
195
196 for (auto it=tCols.begin(); it != tCols.end(); it++)
197 fColMap.insert(*it);
198
199 const fits::Table::Keys& tkeys = fFile->getKeys();
200
201 for (auto it=tkeys.begin(); it != tkeys.end(); it++)
202 fKeyMap.insert(*it);
203
204 return true;
205}
206
207bool FitsDumper::OpenTable(const string &)
208{
209 if (!fFile)
210 {
211 cerr << "No file open." << endl;
212 return false;
213 }
214
215
216 return true;
217}
218
219
220void FitsDumper::List()
221{
222 if (!fFile)
223 {
224 cerr << "No file open." << endl;
225 return;
226 }
227
228 cout << "\nFile: " << fFilename << "\n";
229
230 cout << " " << fKeyMap.find("EXTNAME")->second.value << " [";
231 cout << fKeyMap.find("NAXIS2")->second.value << "]\n";
232
233 for (auto it = fColMap.begin(); it != fColMap.end(); it++)
234 {
235 cout << " " << it->first << "[" << it->second.num << "] (" << it->second.unit << ":" << ValueTypeToStr(it->second.type) << ") ";
236 for (auto jt = fKeyMap.begin(); jt != fKeyMap.end(); jt++)
237 if (jt->second.value == it->first)
238 cout << jt->second.comment << endl;
239 }
240
241 cout << endl;
242 cout << flush;
243}
244
245void FitsDumper::ListKeywords(ostream &out)
246{
247 for (auto it=fKeyMap.begin(); it != fKeyMap.end(); it++) {
248 out << "## " << setw(8) << it->first << " = " << setw(10);
249 out << "'" << it->second.value << "'" << " / " << it->second.comment << endl;
250 }
251}
252
253void FitsDumper::ListHeader()
254{
255 if (!fFile)
256 {
257 cerr << "No table open." << endl;
258 return;
259 }
260
261 cout << "\nTable: " << fKeyMap.find("EXTNAME")->second.value << " (rows=" << fKeyMap.find("NAXIS2")->second.value << ")\n";
262 if (fKeyMap.find("COMMENT") != fKeyMap.end())
263 cout << "Comment: \t" << fKeyMap.find("COMMENT")->second.value << "\n";
264
265 ListKeywords(cout);
266 cout << endl;
267
268}
269
270bool FitsDumper::separateColumnsFromRanges(const vector<string>& list,
271 vector<pair<int, int> >& ranges,
272 vector<string>& listNamesOnly)
273{
274 for (vector<string>::const_iterator it=list.begin(); it!=list.end(); it++)
275 {
276 string columnNameOnly = *it;
277 unsigned long bracketIndex0 = columnNameOnly.find_first_of('[');
278 unsigned long bracketIndex1 = columnNameOnly.find_first_of(']');
279 unsigned long colonIndex = columnNameOnly.find_first_of(':');
280 int columnStart = -1;
281 int columnEnd = -1;
282 if (bracketIndex0 != string::npos)
283 {//there is a range given. Extract the range
284 if (colonIndex != string::npos)
285 {//we have a range here
286 columnStart = atoi(columnNameOnly.substr(bracketIndex0+1, colonIndex-(bracketIndex0+1)).c_str());
287 columnEnd = atoi(columnNameOnly.substr(colonIndex+1, bracketIndex1-(colonIndex+1)).c_str());
288 columnEnd++;
289 }
290 else
291 {//only a single index there
292 columnStart = atoi(columnNameOnly.substr(bracketIndex0+1, bracketIndex1 - (bracketIndex0+1)).c_str());
293 columnEnd = columnStart+1;
294 }
295 columnNameOnly = columnNameOnly.substr(0, bracketIndex0);
296 }
297
298 if (fColMap.find(columnNameOnly) == fColMap.end())
299 {
300 cerr << "ERROR - Column '" << columnNameOnly << "' not found in table." << endl;
301 return false;
302 }
303 fits::Table::Column& cCol = fColMap.find(columnNameOnly)->second;
304 if (bracketIndex0 == string::npos)
305 {//no range given: use the full range
306 ranges.push_back(make_pair(0, cCol.num));
307 columnStart = 0;
308 columnEnd = 1;
309 }
310 else
311 {//use the range extracted earlier
312 if (columnStart < 0)
313 {
314 cerr << "ERROR - Start range for column " << columnNameOnly << " is less than zero (" << columnStart << "). Aborting" << endl;
315 return false;
316 }
317 if (columnEnd>1 && columnEnd > (int)(cCol.num))
318 {
319 cerr << "ERROR - End range for column " << columnNameOnly << " is greater than the last element (" << cCol.num << " vs " << columnEnd << "). Aborting" << endl;
320 return false;
321 }
322 ranges.push_back(make_pair(columnStart, columnEnd));
323 }
324// cout << "Will be exporting from " << columnStart << " to " << columnEnd-1 << " for column " << columnNameOnly << endl;
325 listNamesOnly.push_back(columnNameOnly);
326 }
327 return true;
328}
329// --------------------------------------------------------------------------
330//
331//! Perform the actual dump, based on the current parameters
332//
333bool FitsDumper::Dump(const string &filename, const vector<string> &list, int precision)
334{
335
336 //first of all, let's separate the columns from their ranges and check that the requested columns are indeed part of the file
337 vector<pair<int, int> > ranges;
338 vector<string> listNamesOnly;
339
340 if (!separateColumnsFromRanges(list, ranges, listNamesOnly))
341 {
342 cerr << "Something went wrong while extracting the columns names from parameters. Aborting" << endl;
343 return false;
344 }
345
346 ofstream out(filename=="-"?"/dev/stdout":filename);
347 if (!out)
348 {
349 cerr << "Cannot open file " << filename << ": " << strerror(errno) << endl;
350 return false;
351 }
352
353 out.precision(precision);
354
355 out << "## --------------------------------------------------------------------------\n";
356 if (filename!="-")
357 out << "## File: \t" << filename << '\n';
358 out << "## Table: \t" << fKeyMap.find("EXTNAME")->second.value << '\n';
359 out << "## Comment: \t" << ((fKeyMap.find("COMMENT") != fKeyMap.end()) ? fKeyMap.find("COMMENT")->second.value : "") << '\n';
360
361 out << "## NumRows: \t" << fFile->GetInt("NAXIS2") << '\n';
362 out << "## --------------------------------------------------------------------------\n";
363 ListKeywords(out);
364 out << "## --------------------------------------------------------------------------\n";
365 out << "#\n";
366
367 auto rangesIt = ranges.begin();
368 for (vector<string>::const_iterator it=listNamesOnly.begin(); it!=listNamesOnly.end(); it++, rangesIt++)
369 {
370 const fits::Table::Column& col = fColMap[*it];
371// const MyColumn *col = static_cast<MyColumn*>(fTable->column()[*it]);
372 if (rangesIt->first != 0 || rangesIt->second != (int)(col.num))
373 {
374 out << "#";
375 for (int i=rangesIt->first; i<rangesIt->second; i++)
376 out << " " << *it << "[" << i << "]";
377 out << ": " << col.unit;
378 }
379 else
380 out << "# " << *it << "[" << col.num << "]: " << col.unit;
381
382// FIXME: retrive the column comment
383// if (!col->comment().empty())
384// out << " (" <<col->comment() << ")";
385 out << '\n';
386 }
387 out << "#" << endl;
388
389
390 // Loop over all columns in our list of requested columns
391 vector<pair<char, char*> > columnsData;
392
393 for (vector<string>::const_iterator it=listNamesOnly.begin(); it!=listNamesOnly.end(); it++)
394 {
395 fits::Table::Column& cCol = fColMap.find(*it)->second;
396 columnsData.push_back(make_pair(cCol.type, new char[cCol.num*cCol.size]));
397 fFile->SetPtrAddress(*it, columnsData[columnsData.size()-1].second);
398 }
399 int numRows = fFile->GetInt("NAXIS2");
400 int row = 0;
401 while (fFile->GetNextRow() && row < numRows)
402 {
403 row++;
404 rangesIt = ranges.begin();
405 for (auto it=columnsData.begin(); it != columnsData.end(); it++, rangesIt++)
406 {
407 for (int i=rangesIt->first; i<rangesIt->second; i++)
408 {
409 switch (it->first) {
410 case 'L':
411 out << reinterpret_cast<bool*>(it->second)[i] << " ";
412 break;
413 case 'B':
414 out << reinterpret_cast<char*>(it->second)[i] << " ";
415 break;
416 case 'I':
417 out << reinterpret_cast<int16_t*>(it->second)[i] << " ";
418 break;
419 case 'J':
420 out << reinterpret_cast<int32_t*>(it->second)[i] << " ";
421 break;
422 case 'K':
423 out << reinterpret_cast<int64_t*>(it->second)[i] << " ";
424 break;
425 case 'E':
426 out << reinterpret_cast<float*>(it->second)[i] << " ";
427 break;
428 case 'D':
429 out << reinterpret_cast<double*>(it->second)[i] << " ";
430 break;
431 default:
432 ;
433 }
434 }
435 }
436 }
437 for (auto it = columnsData.begin(); it != columnsData.end(); it++)
438 delete[] it->second;
439 return true;
440}
441
442// --------------------------------------------------------------------------
443//
444//! Retrieves the configuration parameters
445//! @param conf
446//! the configuration object
447//
448int FitsDumper::ExecConfig(Configuration& conf)
449{
450 if (conf.Has("fitsfile"))
451 {
452 if (!OpenFile(conf.Get<string>("fitsfile")))
453 return -1;
454 }
455
456 if (conf.Get<bool>("stat") && conf.Get<bool>("graph"))
457 {
458 cout << "Invalid conbination of options: cannot graph stats. Aborting" << endl;
459 return -1;
460 }
461
462 if (conf.Get<bool>("list"))
463 List();
464
465 if (conf.Has("tablename"))
466 {
467 if (!OpenTable(conf.Get<string>("tablename")))
468 return -1;
469 }
470
471 if (conf.Get<bool>("stat"))
472 {
473 if (!conf.Has("col"))
474 {
475 cout << "Please specify the columns that should be dumped as arguments. Aborting" << endl;
476 return 0;
477 }
478 doStatsPlease(conf.Get<string>("outfile"), conf.Get<vector<string>>("col"), conf.Get<int>("precision"));
479 return 0;
480 }
481
482#ifdef PLOTTING_PLEASE
483 if (conf.Get<bool>("graph"))
484 {
485 if (!conf.Has("col"))
486 {
487 cout << "Please specify the columns that should be dumped as arguments. Aborting" << endl;
488 return 0;
489 }
490 if (conf.Get<bool>("dots"))
491 fDotsPlease = true;
492 doCurvesDisplay(conf.Get<vector<string>>("col"),
493 conf.Get<string>("tablename"));
494 return 1;
495 }
496#endif
497
498 if (conf.Get<bool>("header"))
499 ListHeader();
500
501 if (conf.Get<bool>("header") || conf.Get<bool>("list"))
502 return 1;
503
504 if (conf.Has("outfile"))
505 {
506 if (!conf.Has("col"))
507 {
508 cout << "Please specify the columns that should be dumped as arguments. Aborting" << endl;
509 return 0;
510 }
511 if (!Dump(conf.Get<string>("outfile"),
512 conf.Get<vector<string>>("col"),
513 conf.Get<int>("precision")))
514 return -1;
515 }
516
517
518 return 0;
519}
520
521void PrintUsage()
522{
523 cout <<
524 "fitsdump is a tool to dump data from a FITS table as ascii.\n"
525 "\n"
526 "Usage: fitsdump [OPTIONS] fitsfile col col ... \n"
527 " or: fitsdump [OPTIONS]\n";
528 cout << endl;
529}
530
531void PrintHelp()
532{
533 //
534}
535
536template<typename T>
537void displayStats(char* array, int numElems)
538{
539 sort(reinterpret_cast<T*>(array), &reinterpret_cast<T*>(array)[numElems]);
540 cout << "Min: " << reinterpret_cast<T*>(array)[0] << endl;
541 cout << "Max: " << reinterpret_cast<T*>(array)[numElems-1] << endl;
542 if (numElems%2 == 0)
543 cout << "Med: " << (reinterpret_cast<T*>(array)[numElems/2] + reinterpret_cast<T*>(array)[numElems/2+1])/2.f << endl;
544 else
545 cout << "Med: " << reinterpret_cast<T*>(array)[numElems/2+1] << endl;
546
547 double average = 0;
548 for (int i=0;i<numElems;i++)
549 average += reinterpret_cast<T*>(array)[i];
550 cout << "Mea: " << average/(double)numElems << endl;
551
552}
553int FitsDumper::doStatsPlease(const string &filename, const vector<string>& list, int precision)
554{
555
556 //first of all, let's separate the columns from their ranges and check that the requested columns are indeed part of the file
557 vector<pair<int, int> > ranges;
558 vector<string> listNamesOnly;
559
560 if (!separateColumnsFromRanges(list, ranges, listNamesOnly))
561 {
562 cerr << "Something went wrong while extracting the columns names from parameters. Aborting" << endl;
563 return false;
564 }
565
566 ofstream out(filename=="-"?"/dev/stdout":filename);
567 if (!out)
568 {
569 cerr << "Cannot open file " << filename << ": " << strerror(errno) << endl;
570 return false;
571 }
572
573 out.precision(precision);
574
575 // Loop over all columns in our list of requested columns
576 vector<pair<char, char*> > columnsData;
577 vector<char*> statData;
578 int numRows = fFile->GetInt("NAXIS2");
579 auto rangesIt = ranges.begin();
580 for (vector<string>::const_iterator it=listNamesOnly.begin(); it!=listNamesOnly.end(); it++, rangesIt++)
581 {
582 fits::Table::Column& cCol = fColMap.find(*it)->second;
583 columnsData.push_back(make_pair(cCol.type, new char[cCol.num*cCol.size]));
584 statData.push_back(new char[(rangesIt->second - rangesIt->first)*cCol.size*numRows]);
585 fFile->SetPtrAddress(*it, columnsData[columnsData.size()-1].second);
586 }
587
588 int row = 0;
589 while (fFile->GetNextRow() && row < numRows)
590 {
591 rangesIt = ranges.begin();
592 auto statsIt = statData.begin();
593 for (auto it=columnsData.begin(); it != columnsData.end(); it++, rangesIt++, statsIt++)
594 {
595 int span = rangesIt->second - rangesIt->first;
596 for (int i=rangesIt->first; i<rangesIt->second; i++)
597 {
598 switch (it->first) {
599 case 'L':
600 reinterpret_cast<bool*>(*statsIt)[i - rangesIt->first + row*span] = reinterpret_cast<bool*>(it->second)[i];
601 break;
602 case 'B':
603 reinterpret_cast<char*>(*statsIt)[i - rangesIt->first + row*span] = reinterpret_cast<char*>(it->second)[i];
604 break;
605 case 'I':
606 reinterpret_cast<int16_t*>(*statsIt)[i - rangesIt->first + row*span] = reinterpret_cast<int16_t*>(it->second)[i];
607 break;
608 case 'J':
609 reinterpret_cast<int32_t*>(*statsIt)[i - rangesIt->first + row*span] = reinterpret_cast<int32_t*>(it->second)[i];
610 break;
611 case 'K':
612 reinterpret_cast<int64_t*>(*statsIt)[i - rangesIt->first + row*span] = reinterpret_cast<int64_t*>(it->second)[i];
613 break;
614 case 'E':
615 reinterpret_cast<float*>(*statsIt)[i - rangesIt->first + row*span] = reinterpret_cast<float*>(it->second)[i];
616 break;
617 case 'D':
618 reinterpret_cast<double*>(*statsIt)[i - rangesIt->first + row*span] = reinterpret_cast<double*>(it->second)[i];
619 break;
620 default:
621 ;
622 }
623 }
624 }
625 row++;
626 }
627 for (auto it = columnsData.begin(); it != columnsData.end(); it++)
628 delete[] it->second;
629
630 //okay. So now I've got ALL the data, loaded.
631 //let's do the summing and averaging in a safe way (i.e. avoid overflow of variables as much as possible)
632 rangesIt = ranges.begin();
633 auto statsIt = statData.begin();
634
635 auto nameIt = listNamesOnly.begin();
636 for (auto it=columnsData.begin(); it != columnsData.end(); it++, rangesIt++, statsIt++, nameIt++)
637 {
638 int span = rangesIt->second - rangesIt->first;
639 cout << *nameIt << ": " << endl;
640 switch (it->first) {
641 case 'L':
642 displayStats<bool>(*statsIt, numRows*span);
643 break;
644 case 'B':
645 displayStats<char>(*statsIt, numRows*span);
646 break;
647 case 'I':
648 displayStats<int16_t>(*statsIt, numRows*span);
649 break;
650 case 'J':
651 displayStats<int32_t>(*statsIt, numRows*span);
652 break;
653 case 'K':
654 displayStats<int64_t>(*statsIt, numRows*span);
655 break;
656 case 'E':
657 displayStats<float>(*statsIt, numRows*span);
658 break;
659 case 'D':
660 displayStats<double>(*statsIt, numRows*span);
661 break;
662 default:
663 ;
664 }
665 }
666 return true;
667}
668#ifdef PLOTTING_PLEASE
669int FitsDumper::doCurvesDisplay( const vector<string> &list, const string& tableName)
670{
671 //first of all, let's separate the columns from their ranges and check that the requested columns are indeed part of the file
672 vector<pair<int, int> > ranges;
673 vector<string> listNamesOnly;
674 if (!separateColumnsFromRanges(list, ranges, listNamesOnly))
675 {
676 cerr << "Something went wrong while extracting the columns names from parameters. Aborting" << endl;
677 return false;
678 }
679 vector<string> curvesNames;
680 stringstream str;
681 for (auto it=ranges.begin(), jt=listNamesOnly.begin(); it != ranges.end(); it++, jt++)
682 {
683 for (int i=it->first; i<it->second;i++)
684 {
685 str.str("");
686 str << *jt << "[" << i << "]";
687 curvesNames.push_back(str.str());
688 }
689 }
690 char* handle = new char[17];
691 sprintf(handle,"FitsDump Display");
692// Qt::HANDLE h = *handle;//NULL
693 int argc = 1;
694 char ** argv = &handle;
695 QApplication a(argc, argv);
696
697
698
699 QwtPlot* plot = new QwtPlot();
700 QwtPlotGrid* grid = new QwtPlotGrid;
701 grid->enableX(false);
702 grid->enableY(true);
703 grid->enableXMin(false);
704 grid->enableYMin(false);
705 grid->setMajPen(QPen(Qt::black, 0, Qt::DotLine));
706 grid->attach(plot);
707 plot->setAutoReplot(true);
708 string title = tableName;
709 plot->setAxisScaleDraw( QwtPlot::xBottom, new TimeScaleDraw());
710
711 QWidget window;
712 QHBoxLayout* layout = new QHBoxLayout(&window);
713 layout->setContentsMargins(0,0,0,0);
714 layout->addWidget(plot);
715
716 QwtPlotZoomer zoom(plot->canvas());
717 zoom.setRubberBandPen(QPen(Qt::gray, 2, Qt::DotLine));
718 zoom.setTrackerPen(QPen(Qt::gray));
719 int totalSize = 0;
720 for (unsigned int i=0;i<list.size();i++)
721 totalSize += ranges[i].second - ranges[i].first;
722
723 vector<QwtPlotCurve*> curves(totalSize);
724 int ii=0;
725 for (auto it = curves.begin(), jt=curvesNames.begin(); it != curves.end(); it++, jt++)
726 {
727 *it = new QwtPlotCurve(jt->c_str());
728 switch (ii%6)
729 {
730 case 0:
731 (*it)->setPen(QColor(255,0,0));
732 break;
733 case 1:
734 (*it)->setPen(QColor(0,255,0));
735 break;
736 case 2:
737 (*it)->setPen(QColor(0,0,255));
738 break;
739 case 3:
740 (*it)->setPen(QColor(255,255,0));
741 break;
742 case 4:
743 (*it)->setPen(QColor(0,255,255));
744 break;
745 case 5:
746 (*it)->setPen(QColor(255,0,255));
747 break;
748 default:
749 (*it)->setPen(QColor(0,0,0));
750 };
751 ii++;
752 if (fDotsPlease)
753 (*it)->setStyle(QwtPlotCurve::Dots);
754 else
755 (*it)->setStyle(QwtPlotCurve::Lines);
756 (*it)->attach(plot);
757 }
758 plot->insertLegend(new QwtLegend(), QwtPlot::RightLegend);
759
760
761 vector<pair<char, char*> > columnsData;
762
763 for (vector<string>::const_iterator it=listNamesOnly.begin(); it!=listNamesOnly.end(); it++)
764 {
765 fits::Table::Column& cCol = fColMap.find(*it)->second;
766 columnsData.push_back(make_pair(cCol.type, new char[cCol.num*cCol.size]));
767 fFile->SetPtrAddress(*it, columnsData[columnsData.size()-1].second);
768 }
769
770 //add the time column to the given columns
771 if (fColMap.find("Time") == fColMap.end() && fColMap.find("UnixTimeUTC") == fColMap.end())
772 {
773 cerr << "Error: time column could not be found in given table. Aborting" << endl;
774 return false;
775 }
776 const fits::Table::Column& timeCol = (fColMap.find("Time") != fColMap.end()) ? fColMap.find("Time")->second : fColMap.find("UnixTimeUTC")->second;
777 bool unixTime = (fColMap.find("Time") == fColMap.end());
778 if (unixTime)
779 ranges.push_back(make_pair(0,2));
780 else
781 ranges.push_back(make_pair(0,1));
782 columnsData.push_back(make_pair(timeCol.type, new char[timeCol.num*timeCol.size]));
783 fFile->SetPtrAddress(unixTime ? "UnixTimeUTC" : "Time", columnsData[columnsData.size()-1].second);
784
785// stringstream str;
786 str.str("");
787
788
789 vector<double*> xValues(totalSize);
790 double* yValues;
791 cout.precision(10);
792 str.precision(20);
793 for (auto it=xValues.begin(); it!=xValues.end(); it++)
794 *it = new double[fFile->GetInt("NAXIS2")];
795
796 yValues = new double[fFile->GetInt("NAXIS2")];
797
798 cout.precision(3);
799 int endIndex = 0;
800 int numRows = fFile->GetInt("NAXIS2");
801 for (int i=1;i<numRows;i++)
802 {
803 fFile->GetNextRow();
804 cout << "\r" << "Constructing graph " << ((float)(endIndex)/(float)(fFile->GetInt("NAXIS2")))*100.0 << "%";
805 endIndex++;
806 auto rangesIt = ranges.begin();
807 for (auto it=columnsData.begin(); it != columnsData.end(); it++, rangesIt++)
808 {
809 for (int j=rangesIt->first; j<rangesIt->second; j++)
810 {
811 switch (it->first) {
812 case 'L':
813 str << reinterpret_cast<bool*>(it->second)[j] << " ";
814 break;
815 case 'B':
816 str << reinterpret_cast<char*>(it->second)[j] << " ";
817 break;
818 case 'I':
819 str << reinterpret_cast<int16_t*>(it->second)[j] << " ";
820 break;
821 case 'J':
822 str << reinterpret_cast<int32_t*>(it->second)[j] << " ";
823 break;
824 case 'K':
825 str << reinterpret_cast<int64_t*>(it->second)[j] << " ";
826 break;
827 case 'E':
828 str << reinterpret_cast<float*>(it->second)[j] << " ";
829 break;
830 case 'D':
831 str << reinterpret_cast<double*>(it->second)[j] << " ";
832 break;
833 default:
834 ;
835 }
836 }
837 }
838
839 for (auto it=xValues.begin(); it!= xValues.end(); it++)
840 {
841 str >> (*it)[i-1];
842 }
843 if (unixTime)
844 {
845 long u1, u2;
846 str >> u1 >> u2;
847
848 boost::posix_time::ptime unixTimeT( boost::gregorian::date(1970, boost::gregorian::Jan, 1),
849 boost::posix_time::seconds(u1) + boost::posix_time::microsec(u2));
850
851 Time mjdTime(unixTimeT);
852 yValues[i-1] = mjdTime.Mjd();
853 }
854 else
855 {
856 str >> yValues[i-1];
857 if (yValues[i-1] < 40587)
858 yValues[i-1] += 40587;
859
860 Time t(yValues[i-1]);
861 string time = t.GetAsStr("%H:%M:%S%F");
862 while (time[time.size()-1] == '0' && time.size() > 2)
863 {
864 time = time.substr(0, time.size()-1);
865 }
866 }
867 if (i==1)
868 {
869 Time t(yValues[0]);
870 title += " - " + t.GetAsStr("%Y-%m-%d");
871 plot->setTitle(title.c_str());
872 }
873 }
874 //set the actual data.
875 auto jt = xValues.begin();
876 for (auto it=curves.begin(); it != curves.end(); it++, jt++)
877 (*it)->setRawData(yValues, *jt, endIndex-1);
878
879 QStack<QRectF> stack;
880 double minX, minY, maxX, maxY;
881 minX = minY = 1e10;
882 maxX = maxY = -1e10;
883 QRectF rect;
884 QPointF point;
885 for (auto it=curves.begin(); it!= curves.end(); it++)
886 {
887 rect = (*it)->boundingRect();
888 point = rect.bottomRight();
889 if (point.x() < minX) minX = point.x();
890 if (point.y() < minY) minY = point.y();
891 if (point.x() > maxX) maxX = point.x();
892 if (point.y() > maxY) maxY = point.y();
893 point = rect.topLeft();
894 if (point.x() < minX) minX = point.x();
895 if (point.y() < minY) minY = point.y();
896 if (point.x() > maxX) maxX = point.x();
897 if (point.y() > maxY) maxY = point.y();
898 }
899 QPointF bottomRight(maxX, minY);
900 QPointF topLeft(minX, maxY);
901 QPointF center((bottomRight+topLeft)/2.f);
902 stack.push(QRectF(topLeft + (topLeft-center)*(.5f),bottomRight + (bottomRight-center)*(.5f)));
903 zoom.setZoomStack(stack);
904
905// delete[] fitsBuffer;
906 for (auto it = columnsData.begin(); it != columnsData.end(); it++)
907 delete[] it->second;
908 window.resize(600, 400);
909 window.show();
910
911 a.exec();
912
913
914 for (auto it = curves.begin(); it != curves.end(); it++)
915 {
916 (*it)->detach();
917 delete *it;
918 }
919 grid->detach();
920 for (auto it = xValues.begin(); it != xValues.end(); it++)
921 delete[] *it;
922
923 delete[] yValues;
924
925 delete[] handle;
926
927 return 0;
928}
929#endif
930
931void SetupConfiguration(Configuration& conf)
932{
933 po::options_description configs("Fitsdump options");
934 configs.add_options()
935 ("fitsfile,f", var<string>()
936#if BOOST_VERSION >= 104200
937 ->required()
938#endif
939 , "Name of FITS file")
940 ("tablename,t", var<string>("DATA")
941#if BOOST_VERSION >= 104200
942 ->required()
943#endif
944 , "Name of input table")
945 ("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")
946 ("outfile,o", var<string>("/dev/stdout"), "Name of output file (-:/dev/stdout)")
947 ("precision,p", var<int>(20), "Precision of ofstream")
948 ("list,l", po_switch(), "List all tables and columns in file")
949 ("header,h", po_switch(), "Dump header of given table")
950 ("stat,s", po_switch(), "Perform statistics instead of dump")
951#ifdef PLOTTING_PLEASE
952 ("graph,g", po_switch(), "Plot the columns instead of dumping them")
953 ("dots,d", po_switch(), "Plot using dots instead of lines")
954#endif
955 ;
956
957 po::positional_options_description p;
958 p.add("fitsfile", 1); // The first positional options
959 p.add("col", -1); // All others
960
961 conf.AddOptions(configs);
962 conf.SetArgumentPositions(p);
963}
964
965int main(int argc, const char** argv)
966{
967 Configuration conf(argv[0]);
968 conf.SetPrintUsage(PrintUsage);
969 SetupConfiguration(conf);
970
971 if (!conf.DoParse(argc, argv, PrintHelp))
972 return -1;
973
974 FitsDumper loader;
975 return loader.ExecConfig(conf);
976}
Note: See TracBrowser for help on using the repository browser.