Context Navigation

source: trunk/Mars/mcore/ofits.h@ 17752

Visit:

Last change on this file since 17752 was 17751, checked in by tbretz, 11 years ago
I realized that we were one level too high with deriving the class. Not only there was no possibility to use open/close, furthermore the close funtion was not virtual and therefore not necessarily correctly called. By writing a ofstream type class, both is now solved.
File size: 34.3 KB

Line
1	#ifndef MARS_ofits
2	#define MARS_ofits
3
4	#include <string>
5	#include <string.h>
6	#include <algorithm>
7	#include <sstream>
8	#include <iostream>
9	#include <fstream>
10	#include <iomanip>
11	#include <vector>
12	#include <algorithm>
13	#include <stdexcept>
14
15	#define GCC_VERSION (__GNUC__ * 10000 + __GNUC_MINOR__ * 100 + __GNUC_PATCHLEVEL__)
16
17	#ifdef __CINT__
18	#define off_t size_t
19	#endif
20
21	#ifndef __MARS__
22	#ifndef gLog
23	#define gLog std::cerr
24	#define ___err___ ""
25	#define ___warn___ ""
26	#define ___all___ ""
27	#endif
28	#else
29	#include "MLog.h"
30	#include "MLogManip.h"
31	#define ___err___ err
32	#define ___warn___ warn
33	#define ___all___ all
34	#endif
35
36	#include "FITS.h"
37	#include "checksum.h"
38
39	// Sloppy: allow / <--- left
40	// allow all characters (see specs for what is possible)
41
42	// units: m kg s rad sr K A mol cd Hz J W V N Pa C Ohm S F Wb T Hlm lx
43
44	class ofits : public std::ostream
45	{
46	protected:
47	std::filebuf fFilebuf;
48
49	public:
50	struct Key
51	{
52	std::string key;
53	bool delim;
54	std::string value;
55	std::string comment;
56	std::string fitsString;
57
58	off_t offset; // offset in file
59
60	bool changed; // For closing the file
61
62	Key(const std::string &k="") : key(k), delim(false), fitsString(""), offset(0), changed(true) { }
63
64	std::string Trim(const std::string &str)
65	{
66	// Trim Both leading and trailing spaces
67	const size_t first = str.find_first_not_of(' '); // Find the first character position after excluding leading blank spaces
68	const size_t last = str.find_last_not_of(' '); // Find the first character position from reverse af
69
70	// if all spaces or empty return an empty string
71	if (std::string::npos==first \|\| std::string::npos==last)
72	return std::string();
73
74	return str.substr(first, last-first+1);
75	}
76
77	bool FormatKey()
78	{
79	key = Trim(key);
80	if (key.empty())
81	{
82	#ifdef __EXCEPTIONS
83	throw std::runtime_error("Key name empty.");
84	#else
85	gLog << ___err___ << "ERROR - Key name empty." << std::endl;
86	return false;
87	#endif
88	}
89	if (key.size()>8)
90	{
91	std::ostringstream sout;
92	sout << "Key '" << key << "' exceeds 8 bytes.";
93	#ifdef __EXCEPTIONS
94	throw std::runtime_error(sout.str());
95	#else
96	gLog << ___err___ << "ERROR - " << sout.str() << std::endl;
97	return false;
98	#endif
99	}
100
101	//transform(key.begin(), key.end(), key.begin(), toupper);
102	#if GCC_VERSION < 40603
103	for (std::string::const_iterator c=key.begin(); c<key.end(); c++)
104	#else
105	for (std::string::const_iterator c=key.cbegin(); c<key.cend(); c++)
106	#endif
107	if ((c<'A' \|\| c>'Z') && (c<'0' \|\| c>'9') && c!='-' && c!='_')
108	{
109	std::ostringstream sout;
110	sout << "Invalid character '" << *c << "' found in key '" << key << "'";
111	#ifdef __EXCEPTIONS
112	throw std::runtime_error(sout.str());
113	#else
114	gLog << ___err___ << "ERROR - " << sout.str() << std::endl;
115	return false;
116	#endif
117	}
118
119	return true;
120	}
121
122	bool FormatComment()
123	{
124	comment = Trim(comment);
125
126	#if GCC_VERSION < 40603
127	for (std::string::const_iterator c=key.begin(); c<key.end(); c++)
128	#else
129	for (std::string::const_iterator c=key.cbegin(); c<key.cend(); c++)
130	#endif
131	if (c<32 \|\| c>126)
132	{
133	std::ostringstream sout;
134	sout << "Invalid character '" << c << "' [" << int(c) << "] found in comment '" << comment << "'";
135	#ifdef __EXCEPTIONS
136	throw std::runtime_error(sout.str());
137	#else
138	gLog << ___err___ << "ERROR - " << sout.str() << std::endl;
139	return false;
140	#endif
141	}
142
143	return true;
144	}
145
146	bool check(bool trim=false)
147	{
148	if (!FormatKey())
149	return false;
150
151	if (!FormatComment())
152	return false;
153
154	size_t sz = CalcSize();
155	if (sz<=80)
156	return true;
157
158	if (!trim)
159	{
160	std::ostringstream sout;
161	sout << "Size " << sz << " of entry for key '" << key << "' exceeds 80 characters.";
162	#ifdef __EXCEPTIONS
163	throw std::runtime_error(sout.str());
164	#else
165	gLog << ___err___ << "ERROR - " << sout.str() << std::endl;
166	#endif
167	return false;
168	}
169
170	//looks like something went wrong. Maybe entry is too long ?
171	//try to remove the comment
172	comment = "";
173
174	sz = CalcSize();
175	if (sz<=80)
176	{
177	#ifndef __EXCEPTIONS
178	std::ostringstream sout;
179	sout << "Size " << sz << " of entry for key '" << key << "' exceeds 80 characters... removed comment.";
180	gLog << ___warn___ << "WARNING - " << sout.str() << std::endl;
181	#endif
182	return true;
183	}
184
185	std::ostringstream sout;
186	sout << "Size " << sz << " of entry for key '" << key << "' exceeds 80 characters even without comment.";
187	#ifdef __EXCEPTIONS
188	throw std::runtime_error(sout.str());
189	#else
190	gLog << ___err___ << "ERROR - " << sout.str() << std::endl;
191	return false;
192	#endif
193	}
194
195	size_t CalcSize() const
196	{
197	if (!delim)
198	return 10+comment.size();
199
200	return 10 + (value.size()<20?20:value.size()) + 3 + comment.size();
201	}
202
203	std::string Compile()
204	{
205
206	if (fitsString != "")
207	return fitsString;
208
209	std::ostringstream sout;
210	sout << std::left << std::setw(8) << key;
211
212	if (!delim)
213	{
214	sout << " " << comment;
215	return sout.str();
216	}
217
218	sout << "= ";
219	sout << (!value.empty() && value[0]=='\''?std::left:std::right);
220	sout << std::setw(20) << value << std::left;
221
222	if (!comment.empty())
223	sout << " / " << comment;
224
225	return sout.str();
226	}
227
228	Checksum checksum;
229
230	void Out(std::ostream &fout)
231	{
232	if (!changed)
233	return;
234
235	std::string str = Compile();
236	str.insert(str.end(), 80-str.size(), ' ');
237
238	if (offset==0)
239	offset = fout.tellp();
240
241	//cout << "Write[" << offset << "]: " << key << "/" << value << endl;
242
243	fout.seekp(offset);
244	fout << str;
245
246	checksum.reset();
247	checksum.add(str.c_str(), 80);
248
249	changed = false;
250	}
251	/*
252	void Out(ostream &out)
253	{
254	std::string str = Compile();
255
256	str.insert(str.end(), 80-str.size(), ' ');
257
258	out << str;
259	changed = false;
260	}*/
261	};
262
263	private:
264	std::vector<Key> fKeys;
265
266	std::vector<Key>::iterator findkey(const std::string &key)
267	{
268	for (auto it=fKeys.begin(); it!=fKeys.end(); it++)
269	if (key==it->key)
270	return it;
271
272	return fKeys.end();
273	}
274
275	bool Set(const std::string &key="", bool delim=false, const std::string &value="", const std::string &comment="")
276	{
277	// If no delimit add the row no matter if it alread exists
278	if (delim)
279	{
280	// if the row already exists: update it
281	auto it = findkey(key);
282	if (it!=fKeys.end())
283	{
284	it->value = value;
285	it->changed = true;
286	return true;
287	}
288	}
289
290	if (fTable.num_rows>0)
291	{
292	std::ostringstream sout;
293	sout << "No new header key can be defined, rows were already written to the file... ignoring new key '" << key << "'";
294	#ifdef __EXCEPTIONS
295	throw std::runtime_error(sout.str());
296	#else
297	gLog << ___err___ << "ERROR - " << sout.str() << std::endl;
298	return false;
299	#endif
300	}
301
302	Key entry;
303
304	entry.key = key;
305	entry.delim = delim;
306	entry.value = value;
307	entry.comment = comment;
308	entry.offset = 0;
309	entry.changed = true;
310
311	if (!entry.check(fCommentTrimming))
312	return false;
313
314	fKeys.push_back(entry);
315	return true;
316	}
317
318	protected:
319	struct Table
320	{
321	off_t offset;
322
323	size_t bytes_per_row;
324	size_t num_rows;
325	size_t num_cols;
326
327	struct Column
328	{
329	std::string name;
330	size_t offset;
331	size_t num;
332	size_t size;
333	char type;
334	};
335
336	std::vector<Column> cols;
337
338	Table() : offset(0), bytes_per_row(0), num_rows(0), num_cols(0)
339	{
340	}
341	};
342
343
344	Table fTable;
345
346	std::vector<char> fOutputBuffer;
347
348	std::vector<Table::Column>::const_iterator findcol(const std::string &name)
349	{
350	for (auto it=fTable.cols.cbegin(); it!=fTable.cols.cend(); it++)
351	if (name==it->name)
352	return it;
353
354	return fTable.cols.cend();
355	}
356
357	Checksum fDataSum;
358	Checksum fHeaderSum;
359
360	bool fCommentTrimming;
361	bool fManualExtName;
362
363	public:
364	ofits()
365	: std::ostream(), fFilebuf(), fCommentTrimming(false), fManualExtName(false)
366	{
367	init(&fFilebuf);
368	}
369
370	ofits(const char *fname)
371	: std::ostream(), fFilebuf(), fCommentTrimming(false), fManualExtName(false)
372	{
373	init(&fFilebuf);
374	open(fname);
375	}
376
377	#ifdef __GXX_EXPERIMENTAL_CXX0X__
378	ofits(const std::string &fname)
379	: std::ostream(), fCommentTrimming(false), fManualExtName(false)
380	{
381	init(&fFilebuf);
382	open(fname);
383	}
384	#endif
385
386	virtual ~ofits()
387	{
388	if (is_open())
389	close();
390	}
391	/*
392	filebuf *rdbuf() const
393	{
394	return const_cast<filebuf*>(&fFilebuf);
395	}
396	*/
397	bool is_open()
398	{
399	return fFilebuf.is_open();
400	}
401
402	bool is_open() const
403	{
404	return fFilebuf.is_open();
405	}
406
407	virtual void open(const char *filename, bool addEXTNAMEKey=true)
408	{
409	fDataSum = 0;
410	fHeaderSum = 0;
411
412	fTable = Table();
413	fKeys.clear();
414
415	SetStr("XTENSION", "BINTABLE", "binary table extension");
416	SetInt("BITPIX", 8, "8-bit bytes");
417	SetInt("NAXIS", 2, "2-dimensional binary table");
418	SetInt("NAXIS1", 0, "width of table in bytes");
419	SetInt("NAXIS2", 0, "number of rows in table");
420	SetInt("PCOUNT", 0, "size of special data area");
421	SetInt("GCOUNT", 1, "one data group (required keyword)");
422	SetInt("TFIELDS", 0, "number of fields in each row");
423	if (addEXTNAMEKey)
424	SetStr("EXTNAME", "", "name of extension table");
425	else
426	fManualExtName = true;
427	SetStr("CHECKSUM", "0000000000000000", "Checksum for the whole HDU");
428	SetStr("DATASUM", " 0", "Checksum for the data block");
429
430	if (!fFilebuf.open(filename, ios_base::out\|ios_base::trunc))
431	setstate(ios_base::failbit);
432	else
433	clear();
434	}
435
436	virtual void open(const std::string &filename, bool addEXTNAMEKey=true)
437	{
438	open(filename.c_str(), addEXTNAMEKey);
439	}
440
441	void AllowCommentsTrimming(bool allow)
442	{
443	fCommentTrimming = allow;
444	}
445	//Etienne: required to enable 1 to 1 reconstruction of files
446	bool SetKeyComment(const std::string& key, const std::string& comment)
447	{
448	auto it = findkey(key);
449	if (it==fKeys.end())
450	return false;
451	it->comment = comment;
452	it->changed = true;
453	return true;
454	}
455	bool SetKeyFromFitsString(const std::string& fitsString)
456	{
457	if (fTable.num_rows>0)
458	{
459	std::ostringstream sout;
460	sout << "No new header key can be defined, rows were already written to the file... ignoring new key '" << fitsString << "'";
461	#ifdef __EXCEPTIONS
462	throw std::runtime_error(sout.str());
463	#else
464	gLog << ___err___ << "ERROR - " << sout.str() << std::endl;
465	return false;
466	#endif
467	}
468
469	Key entry;
470	entry.fitsString = fitsString;
471	entry.changed = true;
472	fKeys.push_back(entry);
473	return true;
474	}
475	bool SetRaw(const std::string &key, const std::string &val, const std::string &comment)
476	{
477	return Set(key, true, val, comment);
478	}
479
480	bool SetBool(const std::string &key, bool b, const std::string &comment="")
481	{
482	return Set(key, true, b?"T":"F", comment);
483	}
484
485	bool AddEmpty(const std::string &key, const std::string &comment="")
486	{
487	return Set(key, true, "", comment);
488	}
489
490	bool SetStr(const std::string &key, std::string s, const std::string &comment="")
491	{
492	for (uint i=0; i<s.length(); i++)
493	if (s[i]=='\'')
494	s.insert(i++, "\'");
495
496	return Set(key, true, "'"+s+"'", comment);
497	}
498
499	bool SetInt(const std::string &key, int64_t i, const std::string &comment="")
500	{
501	std::ostringstream sout;
502	sout << i;
503
504	return Set(key, true, sout.str(), comment);
505	}
506
507	bool SetFloat(const std::string &key, double f, int p, const std::string &comment="")
508	{
509	std::ostringstream sout;
510
511	if (p<0)
512	sout << std::setprecision(-p) << fixed;
513	if (p>0)
514	sout << std::setprecision(p);
515	if (p==0)
516	sout << std::setprecision(f>1e-100 && f<1e100 ? 15 : 14);
517
518	sout << f;
519
520	std::string str = sout.str();
521
522	replace(str.begin(), str.end(), 'e', 'E');
523
524	if (str.find_first_of('E')==std::string::npos && str.find_first_of('.')==std::string::npos)
525	str += ".";
526
527	return Set(key, true, str, comment);
528	}
529
530	bool SetFloat(const std::string &key, double f, const std::string &comment="")
531	{
532	return SetFloat(key, f, 0, comment);
533	}
534
535	bool SetHex(const std::string &key, uint64_t i, const std::string &comment="")
536	{
537	std::ostringstream sout;
538	sout << std::hex << "0x" << i;
539	return SetStr(key, sout.str(), comment);
540	}
541
542	bool AddComment(const std::string &comment)
543	{
544	return Set("COMMENT", false, "", comment);
545	}
546
547	bool AddHistory(const std::string &comment)
548	{
549	return Set("HISTORY", false, "", comment);
550	}
551
552	void End()
553	{
554	Set("END");
555	while (fKeys.size()%36!=0)
556	fKeys.emplace_back();
557	}
558
559	std::string CommentFromType(char type)
560	{
561	std::string comment;
562
563	switch (type)
564	{
565	case 'L': comment = "[1-byte BOOL]"; break;
566	case 'A': comment = "[1-byte CHAR]"; break;
567	case 'B': comment = "[1-byte BOOL]"; break;
568	case 'I': comment = "[2-byte INT]"; break;
569	case 'J': comment = "[4-byte INT]"; break;
570	case 'K': comment = "[8-byte INT]"; break;
571	case 'E': comment = "[4-byte FLOAT]"; break;
572	case 'D': comment = "[8-byte FLOAT]"; break;
573	case 'Q': comment = "[var. Length]"; break;
574	}
575
576	return comment;
577	}
578
579	uint32_t SizeFromType(char type)
580	{
581	size_t size = 0;
582
583	switch (type)
584	{
585	case 'L': size = 1; break;
586	case 'A': size = 1; break;
587	case 'B': size = 1; break;
588	case 'I': size = 2; break;
589	case 'J': size = 4; break;
590	case 'K': size = 8; break;
591	case 'E': size = 4; break;
592	case 'D': size = 8; break;
593	case 'Q': size = 16; break;
594	}
595
596	return size;
597	}
598	//ETIENNE to be able to restore the file 1 to 1, I must restore the header keys myself
599	virtual bool AddColumn(uint32_t cnt, char typechar, const std::string &name, const std::string &unit, const std::string &comment="", bool addHeaderKeys=true)
600	{
601	if (tellp()<0)
602	{
603	std::ostringstream sout;
604	sout << "File not open... ignoring column '" << name << "'";
605	#ifdef __EXCEPTIONS
606	throw std::runtime_error(sout.str());
607	#else
608	gLog << ___err___ << "ERROR - " << sout.str() << std::endl;
609	return false;
610	#endif
611	}
612
613	if (tellp()>0)
614	{
615	std::ostringstream sout;
616	sout << "Header already written, no new column can be defined... ignoring column '" << name << "'";
617	#ifdef __EXCEPTIONS
618	throw std::runtime_error(sout.str());
619	#else
620	gLog << ___err___ << "ERROR - " << sout.str() << std::endl;
621	return false;
622	#endif
623	}
624
625	if (findcol(name)!=fTable.cols.cend())
626	{
627	std::ostringstream sout;
628	sout << "A column with the name '" << name << "' already exists.";
629	#ifdef __EXCEPTIONS
630	throw std::runtime_error(sout.str());
631	#else
632	gLog << ___err___ << "ERROR - " << sout.str() << std::endl;
633	return false;
634	#endif
635	}
636
637	typechar = toupper(typechar);
638
639	static const std::string allow("LABIJKEDQ");
640	#if GCC_VERSION < 40603
641	if (std::find(allow.begin(), allow.end(), typechar)==allow.end())
642	#else
643	if (std::find(allow.cbegin(), allow.cend(), typechar)==allow.end())
644	#endif
645	{
646	std::ostringstream sout;
647	sout << "Column type '" << typechar << "' not supported.";
648	#ifdef __EXCEPTIONS
649	throw std::runtime_error(sout.str());
650	#else
651	gLog << ___err___ << "ERROR - " << sout.str() << std::endl;
652	return false;
653	#endif
654	}
655
656	std::ostringstream type;
657	type << cnt;
658	if (typechar=='Q')
659	type << "QB";
660	else
661	type << typechar;
662
663	fTable.num_cols++;
664
665	if (addHeaderKeys)
666	{
667	#if GCC_VERSION < 40603
668	const std::string nc = std::to_string((long long int)(fTable.num_cols));
669	#else
670	const std::string nc = std::to_string(fTable.num_cols);
671	#endif
672	SetStr("TFORM"+nc, type.str(), "format of "+name+" "+CommentFromType(typechar));
673	SetStr("TTYPE"+nc, name, comment);
674	if (!unit.empty())
675	SetStr("TUNIT"+nc, unit, "unit of "+name);
676	}
677	size_t size = SizeFromType(typechar);
678
679	Table::Column col;
680
681	col.name = name;
682	col.type = typechar;
683	col.num = cnt;
684	col.size = size;
685	col.offset = fTable.bytes_per_row;
686
687	fTable.cols.push_back(col);
688
689	fTable.bytes_per_row += size*cnt;
690
691	// Align to four bytes
692	fOutputBuffer.resize(fTable.bytes_per_row + 4-fTable.bytes_per_row%4);
693
694	return true;
695	}
696
697	virtual bool AddColumn(const FITS::Compression&, uint32_t cnt, char typechar, const std::string& name, const std::string& unit, const std::string& comment="", bool addHeaderKeys=true)
698	{
699	return AddColumn(cnt, typechar, name, unit, comment, addHeaderKeys);
700	}
701
702	bool AddColumnShort(uint32_t cnt, const std::string &name, const std::string &unit="", const std::string &comment="")
703	{ return AddColumn(cnt, 'I', name, unit, comment); }
704	bool AddColumnInt(uint32_t cnt, const std::string &name, const std::string &unit="", const std::string &comment="")
705	{ return AddColumn(cnt, 'J', name, unit, comment); }
706	bool AddColumnLong(uint32_t cnt, const std::string &name, const std::string &unit="", const std::string &comment="")
707	{ return AddColumn(cnt, 'K', name, unit, comment); }
708	bool AddColumnFloat(uint32_t cnt, const std::string &name, const std::string &unit="", const std::string &comment="")
709	{ return AddColumn(cnt, 'E', name, unit, comment); }
710	bool AddColumnDouble(uint32_t cnt, const std::string &name, const std::string &unit="", const std::string &comment="")
711	{ return AddColumn(cnt, 'D', name, unit, comment); }
712	bool AddColumnChar(uint32_t cnt, const std::string &name, const std::string &unit="", const std::string &comment="")
713	{ return AddColumn(cnt, 'A', name, unit, comment); }
714	bool AddColumnByte(uint32_t cnt, const std::string &name, const std::string &unit="", const std::string &comment="")
715	{ return AddColumn(cnt, 'B', name, unit, comment); }
716	bool AddColumnBool(uint32_t cnt, const std::string &name, const std::string &unit="", const std::string &comment="")
717	{ return AddColumn(cnt, 'L', name, unit, comment); }
718
719	bool AddColumnShort(const std::string &name, const std::string &unit="", const std::string &comment="")
720	{ return AddColumn(1, 'I', name, unit, comment); }
721	bool AddColumnInt(const std::string &name, const std::string &unit="", const std::string &comment="")
722	{ return AddColumn(1, 'J', name, unit, comment); }
723	bool AddColumnLong(const std::string &name, const std::string &unit="", const std::string &comment="")
724	{ return AddColumn(1, 'K', name, unit, comment); }
725	bool AddColumnFloat(const std::string &name, const std::string &unit="", const std::string &comment="")
726	{ return AddColumn(1, 'E', name, unit, comment); }
727	bool AddColumnDouble(const std::string &name, const std::string &unit="", const std::string &comment="")
728	{ return AddColumn(1, 'D', name, unit, comment); }
729	bool AddColumnChar(const std::string &name, const std::string &unit="", const std::string &comment="")
730	{ return AddColumn(1, 'A', name, unit, comment); }
731	bool AddColumnByte(const std::string &name, const std::string &unit="", const std::string &comment="")
732	{ return AddColumn(1, 'B', name, unit, comment); }
733	bool AddColumnBool(const std::string &name, const std::string &unit="", const std::string &comment="")
734	{ return AddColumn(1, 'L', name, unit, comment); }
735
736	bool AddColumnShort(const FITS::Compression &comp, uint32_t cnt, const std::string &name, const std::string &unit="", const std::string &comment="")
737	{ return AddColumn(comp, cnt, 'I', name, unit, comment); }
738	bool AddColumnInt(const FITS::Compression &comp, uint32_t cnt, const std::string &name, const std::string &unit="", const std::string &comment="")
739	{ return AddColumn(comp, cnt, 'J', name, unit, comment); }
740	bool AddColumnLong(const FITS::Compression &comp, uint32_t cnt, const std::string &name, const std::string &unit="", const std::string &comment="")
741	{ return AddColumn(comp, cnt, 'K', name, unit, comment); }
742	bool AddColumnFloat(const FITS::Compression &comp, uint32_t cnt, const std::string &name, const std::string &unit="", const std::string &comment="")
743	{ return AddColumn(comp, cnt, 'E', name, unit, comment); }
744	bool AddColumnDouble(const FITS::Compression &comp, uint32_t cnt, const std::string &name, const std::string &unit="", const std::string &comment="")
745	{ return AddColumn(comp, cnt, 'D', name, unit, comment); }
746	bool AddColumnChar(const FITS::Compression &comp, uint32_t cnt, const std::string &name, const std::string &unit="", const std::string &comment="")
747	{ return AddColumn(comp, cnt, 'A', name, unit, comment); }
748	bool AddColumnByte(const FITS::Compression &comp, uint32_t cnt, const std::string &name, const std::string &unit="", const std::string &comment="")
749	{ return AddColumn(comp, cnt, 'B', name, unit, comment); }
750	bool AddColumnBool(const FITS::Compression &comp, uint32_t cnt, const std::string &name, const std::string &unit="", const std::string &comment="")
751	{ return AddColumn(comp, cnt, 'L', name, unit, comment); }
752
753	bool AddColumnShort(const FITS::Compression &comp, const std::string &name, const std::string &unit="", const std::string &comment="")
754	{ return AddColumn(comp, 1, 'I', name, unit, comment); }
755	bool AddColumnInt(const FITS::Compression &comp, const std::string &name, const std::string &unit="", const std::string &comment="")
756	{ return AddColumn(comp, 1, 'J', name, unit, comment); }
757	bool AddColumnLong(const FITS::Compression &comp, const std::string &name, const std::string &unit="", const std::string &comment="")
758	{ return AddColumn(comp, 1, 'K', name, unit, comment); }
759	bool AddColumnFloat(const FITS::Compression &comp, const std::string &name, const std::string &unit="", const std::string &comment="")
760	{ return AddColumn(comp, 1, 'E', name, unit, comment); }
761	bool AddColumnDouble(const FITS::Compression &comp, const std::string &name, const std::string &unit="", const std::string &comment="")
762	{ return AddColumn(comp, 1, 'D', name, unit, comment); }
763	bool AddColumnChar(const FITS::Compression &comp, const std::string &name, const std::string &unit="", const std::string &comment="")
764	{ return AddColumn(comp, 1, 'A', name, unit, comment); }
765	bool AddColumnByte(const FITS::Compression &comp, const std::string &name, const std::string &unit="", const std::string &comment="")
766	{ return AddColumn(comp, 1, 'B', name, unit, comment); }
767	bool AddColumnBool(const FITS::Compression &comp, const std::string &name, const std::string &unit="", const std::string &comment="")
768	{ return AddColumn(comp, 1, 'L', name, unit, comment); }
769
770	/*
771	bool AddKey(string key, double d, const std::string &comment)
772	{
773	ostringstream out;
774	out << d;
775
776	std::string s = out.str();
777
778	replace(s.begin(), s.end(), "e", "E");
779
780	return AddKey(key, s, comment);
781	}*/
782
783
784	Checksum WriteHeader(std::ostream &fout)
785	{
786	Checksum sum;
787	uint32_t count=0;
788	for (auto it=fKeys.begin(); it!=fKeys.end(); it++)
789	{
790	it->Out(fout);
791	sum += it->checksum;
792	count++;
793	}
794	fout.flush();
795
796	return sum;
797	}
798
799	Checksum WriteHeader()
800	{
801	return WriteHeader(*this);
802	}
803
804	void FlushHeader()
805	{
806	const off_t pos = tellp();
807	WriteHeader();
808	seekp(pos);
809	}
810
811	Checksum WriteFitsHeader()
812	{
813	ofits h;
814
815	h.SetBool("SIMPLE", true, "file does conform to FITS standard");
816	h.SetInt ("BITPIX", 8, "number of bits per data pixel");
817	h.SetInt ("NAXIS", 0, "number of data axes");
818	h.SetBool("EXTEND", true, "FITS dataset may contain extensions");
819	h.SetStr ("CHECKSUM","0000000000000000", "Checksum for the whole HDU");
820	h.SetStr ("DATASUM", " 0", "Checksum for the data block");
821	h.AddComment("FITS (Flexible Image Transport System) format is defined in 'Astronomy");
822	h.AddComment("and Astrophysics', volume 376, page 359; bibcode: 2001A&A...376..359H");
823	h.End();
824
825	const Checksum sum = h.WriteHeader(*this);
826
827	h.SetStr("CHECKSUM", sum.str());
828
829	const size_t offset = tellp();
830	h.WriteHeader(*this);
831	seekp(offset);
832
833	return sum;
834	}
835
836	virtual bool WriteDrsOffsetsTable ()
837	{
838	return true;
839	}
840
841	virtual bool WriteCatalog()
842	{
843	return true;
844	}
845
846	virtual bool WriteTableHeader(const char *name="DATA")
847	{
848	if (tellp()>0)
849	{
850	#ifdef __EXCEPTIONS
851	throw std::runtime_error("File not empty anymore.");
852	#else
853	gLog << ___err___ << "ERROR - File not empty anymore." << std::endl;
854	return false;
855	#endif
856	}
857
858	fHeaderSum = WriteFitsHeader();
859
860	WriteDrsOffsetsTable();
861
862	if (!fManualExtName)
863	SetStr("EXTNAME", name);
864	SetInt("NAXIS1", fTable.bytes_per_row);
865	SetInt("TFIELDS", fTable.cols.size());
866
867	End();
868
869	WriteHeader();
870
871	WriteCatalog();
872
873	return good();
874	}
875
876	template<size_t N>
877	void revcpy(char dest, const char src, int num)
878	{
879	const char pend = src + numN;
880	for (const char *ptr = src; ptr<pend; ptr+=N, dest+=N)
881	std::reverse_copy(ptr, ptr+N, dest);
882	}
883
884	virtual uint32_t GetBytesPerRow() const { return fTable.bytes_per_row; }
885
886	virtual bool WriteRow(const void *ptr, size_t cnt, bool byte_swap=true)
887	{
888	// FIXME: Make sure that header was already written
889	// or write header now!
890	if (cnt!=fTable.bytes_per_row)
891	{
892	std::ostringstream sout;
893	sout << "WriteRow - Size " << cnt << " does not match expected size " << fTable.bytes_per_row;
894	#ifdef __EXCEPTIONS
895	throw std::runtime_error(sout.str());
896	#else
897	gLog << ___err___ << "ERROR - " << sout.str() << std::endl;
898	return false;
899	#endif
900	}
901
902	// For the checksum we need everything to be correctly aligned
903	const uint8_t offset = fTable.offset%4;
904
905	char *buffer = fOutputBuffer.data() + offset;
906
907	auto ib = fOutputBuffer.begin();
908	auto ie = fOutputBuffer.rbegin();
909	*ib++ = 0;
910	*ib++ = 0;
911	*ib++ = 0;
912	*ib = 0;
913
914	*ie++ = 0;
915	*ie++ = 0;
916	*ie++ = 0;
917	*ie = 0;
918
919	if (!byte_swap)
920	memcpy(buffer, ptr, cnt);
921	else
922	{
923	for (auto it=fTable.cols.cbegin(); it!=fTable.cols.cend(); it++)
924	{
925	const char src = reinterpret_cast<const char>(ptr) + it->offset;
926	char *dest = buffer + it->offset;
927
928	// Let the compiler do some optimization by
929	// knowing the we only have 1, 2, 4 and 8
930	switch (it->size)
931	{
932	case 1: memcpy (dest, src, it->num*it->size); break;
933	case 2: revcpy<2>(dest, src, it->num); break;
934	case 4: revcpy<4>(dest, src, it->num); break;
935	case 8: revcpy<8>(dest, src, it->num); break;
936	}
937	}
938	}
939
940	write(buffer, cnt);
941	fDataSum.add(fOutputBuffer);
942
943	fTable.num_rows++;
944	fTable.offset += cnt;
945	return good();
946	}
947
948	template<typename N>
949	bool WriteRow(const std::vector<N> &vec)
950	{
951	return WriteRow(vec.data(), vec.size()*sizeof(N));
952	}
953
954	// Flushes the number of rows to the header on disk
955	virtual void FlushNumRows()
956	{
957	SetInt("NAXIS2", fTable.num_rows);
958	FlushHeader();
959	}
960
961	size_t GetNumRows() const { return fTable.num_rows; }
962
963	void AlignTo2880Bytes()
964	{
965	if (tellp()%(80*36)>0)
966	{
967	std::vector<char> filler(8036-tellp()%(8036));
968	write(filler.data(), filler.size());
969	}
970	}
971
972	Checksum UpdateHeaderChecksum()
973	{
974	std::ostringstream dataSumStr;
975	dataSumStr << fDataSum.val();
976	SetStr("DATASUM", dataSumStr.str());
977
978	const Checksum sum = WriteHeader();
979
980	//sum += headersum;
981
982	SetStr("CHECKSUM", (sum+fDataSum).str());
983
984	return WriteHeader();
985	}
986	virtual bool close()
987	{
988	if (tellp()<0)
989	return false;
990
991	AlignTo2880Bytes();
992
993	// We don't have to jump back to the end of the file
994	SetInt("NAXIS2", fTable.num_rows);
995
996	const Checksum chk = UpdateHeaderChecksum();
997
998	if (!fFilebuf.close())
999	setstate(ios_base::failbit);
1000
1001	if ((chk+fDataSum).valid())
1002	return true;
1003
1004	std::ostringstream sout;
1005	sout << "Checksum (" << std::hex << chk.val() << ") invalid.";
1006	#ifdef __EXCEPTIONS
1007	throw std::runtime_error(sout.str());
1008	#else
1009	gLog << ___err___ << "ERROR - " << sout.str() << std::endl;
1010	return false;
1011	#endif
1012	}
1013
1014	std::pair<std::string, int> GetChecksumData()
1015	{
1016	std::string datasum;
1017	std::string checksum;
1018	//cannot directly use the Get methods, because they are only in fits.h
1019	for (std::vector<Key>::const_iterator it=fKeys.cbegin(); it!= fKeys.cend(); it++)
1020	{
1021	if (it->key == "CHECKSUM") checksum = it->value;
1022	if (it->key == "DATASUM") datasum = it->value;
1023	}
1024	if (checksum[0] == '\'')
1025	checksum = checksum.substr(1,checksum.size()-2);
1026	if (datasum[0] == '\'')
1027	datasum = datasum.substr(1, datasum.size()-2);
1028	return std::make_pair(checksum, atoi(datasum.c_str()));
1029	}
1030
1031	void SetDefaultKeys()
1032	{
1033	SetStr("TELESCOP", "FACT", "Telescope that acquired this data");
1034	SetStr("CREATOR", typeid(*this).name(), "Class that wrote this file");
1035	SetFloat("EXTREL", 1.0, "Release Number");
1036	SetStr("COMPILED", __DATE__" " __TIME__, "Compile time");
1037	SetStr("ORIGIN", "FACT", "Institution that wrote the file");
1038	SetStr("TIMESYS", "UTC", "Time system");
1039	SetStr("TIMEUNIT", "d", "Time given in days w.r.t. to MJDREF");
1040	SetInt("MJDREF", 40587, "MJD to UNIX time (seconds since 1970/1/1)");
1041	SetStr("PACKAGE", PACKAGE_NAME, "Package name");
1042	SetStr("VERSION", PACKAGE_VERSION, "Package description");
1043	SetStr("REVISION", REVISION, "SVN revision");
1044
1045	const time_t t0 = time(NULL);
1046	const struct tm *tmp1 = gmtime(&t0);
1047
1048	std::string str(19, '\0');
1049	if (tmp1 && strftime(const_cast<char*>(str.data()), 20, "%Y-%m-%dT%H:%M:%S", tmp1))
1050	SetStr("DATE", str, "File creation date");
1051	}
1052	};
1053
1054	#if 0
1055	#include "fits.h"
1056
1057	int main()
1058	{
1059	using namespace std;
1060
1061	ofits h2("delme.fits");
1062
1063	h2.SetInt("KEY1", 1, "comment 1");
1064	h2.AddColumnInt(2, "testcol1", "counts", "My comment");
1065	h2.AddColumnInt("testcol2", "counts", "My comment");
1066	//h2.AddColumnInt("testcol2", "counts", "My comment");
1067	h2.SetInt("KEY2", 2, "comment 2");
1068
1069	/*
1070	AddFloat("X0", 0.000123456, "number of fields in each row");
1071	AddFloat("X1", 0, "number of fields in each row");
1072	AddFloat("X2", 12345, "number of fields in each row");
1073	AddFloat("X3", 123456.67890, "number of fields in each row");
1074	AddFloat("X4", 1234567890123456789.12345678901234567890, "number of fields in each row");
1075	AddFloat("X5", 1234567890.1234567890e20, "number of fields in each row");
1076	AddFloat("X6", 1234567890.1234567890e-20, "number of fields in each row");
1077	AddFloat("XB", 1234567890.1234567890e-111, "number of fields in each row");
1078	AddFloat("X7", 1e-5, "number of fields in each row");
1079	AddFloat("X8", 1e-6, "number of fields in each row");
1080	//AddStr("12345678", "123456789012345678", "12345678901234567890123456789012345678901234567");
1081	*/
1082	// -
1083
1084	h2.WriteTableHeader("TABLE_NAME");
1085
1086	for (int i=0; i<10; i++)
1087	{
1088	int j[3] = { i+10, i10, i100 };
1089	h2.WriteRow(j, 3*sizeof(i));
1090	}
1091
1092	//h2.AddColumnInt("testcol2", "counts", "My comment");
1093	//h2.SetInt("KEY3", 2, "comment 2");
1094	h2.SetInt("KEY2", 2, "comment 2xxx");
1095	h2.SetInt("KEY1", 11);
1096
1097	h2.close();
1098
1099	cout << "---" << endl;
1100
1101	fits f("delme.fits");
1102	if (!f)
1103	throw std::runtime_error("xxx");
1104
1105	cout << "Header is valid: " << f.IsHeaderOk() << endl;
1106
1107	while (f.GetNextRow());
1108
1109	cout << "File is valid: " << f.IsFileOk() << endl;
1110
1111	cout << "---" << endl;
1112
1113	return 0;
1114	}
1115	#endif
1116
1117	#endif

Note: See TracBrowser for help on using the repository browser.

Download in other formats: