Logo ROOT  
Reference Guide
 
Loading...
Searching...
No Matches
TFITS.cxx
Go to the documentation of this file.
1// @(#)root/graf2d:$Id$
2// Author: Claudi Martinez, July 19th 2010
3
4/*************************************************************************
5 * Copyright (C) 1995-2010, Rene Brun and Fons Rademakers. *
6 * All rights reserved. *
7 * *
8 * For the licensing terms see $ROOTSYS/LICENSE. *
9 * For the list of contributors see $ROOTSYS/README/CREDITS. *
10 *************************************************************************/
11
12/// \defgroup fitsio FITS file
13/// \brief Interface to FITS file.
14/// \ingroup Graphics2D
15///
16/// TFITS is an interface that lets you reading Flexible Image Transport System
17/// (FITS) files, which are generally used in astronomy. This file format
18/// was standardized 1981 and today is still widely used among professional
19/// and amateur astronomers. FITS is not only an image file, but also
20/// it can contain spectrums, data tables, histograms, and multidimensional
21/// data. Furthermore, FITS data can be described itself by containing
22/// human-readable information that let us to interpret the data within
23/// the FITS file. For example, a FITS could contain a 3D data cube,
24/// but an additional description would tell us that we must read it, for
25/// example, as a 3-layer image.
26///
27/// TFITS requires CFITSIO library to be installed on your system. It
28/// is currently maintained by NASA/GSFC and can be downloaded from
29/// [NASA/GSFC web site](http://fits.gsfc.nasa.gov), as well as documentation.
30///
31/// Using this interface is easy and straightforward. There is only 1 class
32/// called "TFITSHDU" which has several methods to extract data from a
33/// FITS file, more specifically, from an HDU within the file. An HDU, or
34/// Header Data Unit, is a chunk of data with a header containing several
35/// "keyword = value" tokens. The header describes the structure of data
36/// within the HDU. An HDU can be of two types: an "image HDU" or a "table
37/// HDU". The former can be any kind of multidimensional array of real numbers,
38/// by which the name "image" may be confusing: you can store an image, but
39/// you can also store a N-dimensional data cube. On the other hand, table
40/// HDUs are sets of several rows and columns (a.k.a fields) which contain
41/// generic data, as strings, real or complex numbers and even arrays.
42///
43/// Please have a look to the tutorials ($ROOTSYS/tutorials/fitsio/) to see
44/// some examples. IMPORTANT: to run tutorials it is required that
45/// you change the current working directory of ROOT (CINT) shell to the
46/// tutorials directory. Example:
47/// ~~~ {.cpp}
48/// root [1] gSystem->ChangeDirectory("tutorials/fitsio")
49/// root [1] .x FITS_tutorial1.C
50/// ~~~
51/// LIST OF TODO
52/// - Support for complex values within data tables
53/// - Support for reading arrays from table cells
54/// - Support for grouping
55///
56/// IMPLEMENTATION NOTES:
57///
58/// CFITSIO library uses standard C types ('int', 'long', ...). To avoid
59/// confusion, the same types are used internally by the access methods.
60/// However, class's fields are ROOT-defined types.
61
62/** \class TFITSHDU
63\ingroup fitsio
64
65FITS file interface class
66
67TFITSHDU is a class that allows extracting images and data from FITS files and contains
68several methods to manage them.
69*/
70
71#include "TFITS.h"
72#include "TImage.h"
73#include "TArrayI.h"
74#include "TArrayD.h"
75#include "TH1D.h"
76#include "TH2D.h"
77#include "TH3D.h"
78#include "TVectorD.h"
79#include "TMatrixD.h"
80#include "TObjArray.h"
81#include "TObjString.h"
82#include "TCanvas.h"
83#include "TMath.h"
84#include "strlcpy.h"
85
86#include "fitsio.h"
87#include <cstdlib>
88
89
90////////////////////////////////////////////////////////////////////////////////
91/// Clean path from possible filter and put the result in 'dst'.
92
94{
96
97 Ssiz_t ndx = dst.Index("[", 1, 0, TString::kExact);
98 if (ndx != kNPOS) {
99 dst.Resize(ndx);
100 }
101}
102
103
104////////////////////////////////////////////////////////////////////////////////
105/// TFITSHDU constructor from file path with HDU selection filter.
106/// Please refer to CFITSIO manual for more information about
107/// HDU selection filters.
108///
109/// Examples:
110/// - `TFITSHDU("/path/to/myfile.fits")`: just open the PRIMARY HDU
111/// - `TFITSHDU("/path/to/myfile.fits[1]")`: open HDU #1
112/// - `TFITSHDU("/path/to/myfile.fits[PICS]")`: open HDU called 'PICS'
113/// - `TFITSHDU("/path/to/myfile.fits[ACQ][EXPOSURE > 5]")`: open the (table) HDU called 'ACQ' and
114/// selects the rows that have column 'EXPOSURE'
115/// greater than 5.
116
129
130////////////////////////////////////////////////////////////////////////////////
131/// TFITSHDU constructor from filepath and extension number.
132
134{
137
138 //Add "by extension number" filter
140
141 if (kFALSE == LoadHDU(fFilePath)) {
143 throw -1;
144 }
145}
146
147////////////////////////////////////////////////////////////////////////////////
148/// TFITSHDU constructor from filepath and extension name.
149
151{
154
155 //Add "by extension number" filter
157
158
159 if (kFALSE == LoadHDU(fFilePath)) {
161 throw -1;
162 }
163}
164
165////////////////////////////////////////////////////////////////////////////////
166/// TFITSHDU destructor.
167
172
173////////////////////////////////////////////////////////////////////////////////
174/// Release internal resources.
175
177{
178 if (fRecords) delete [] fRecords;
179
180 if (fType == kImageHDU) {
181 if (fSizes) delete fSizes;
182 if (fPixels) delete fPixels;
183 } else {
184 if (fColumnsInfo) {
185 if (fCells) {
186 for (Int_t i = 0; i < fNColumns; i++) {
187 if (fColumnsInfo[i].fType == kString) {
188 //Deallocate character arrays allocated for kString columns
189 Int_t offset = i * fNRows;
190 for (Int_t row = 0; row < fNRows; row++) {
191 delete [] fCells[offset+row].fString;
192 }
193 } else if (fColumnsInfo[i].fType == kRealArray) {
194 //Deallocate character arrays allocated for kString columns
195 Int_t offset = i * fNRows;
196 for (Int_t row = 0; row < fNRows; row++) {
197 delete [] fCells[offset+row].fRealArray;
198 }
199 } else if (fColumnsInfo[i].fType == kRealVector) {
200 // Deallocate character arrays allocated for variable-length array columns
201 Int_t offset = i * fNRows;
202 for (Int_t row = 0; row < fNRows; row++) {
203 delete fCells[offset + row].fRealVector;
204 }
205 }
206 }
207
208 delete [] fCells;
209 }
210
211 delete [] fColumnsInfo;
212 }
213
214
215 }
216}
217
218////////////////////////////////////////////////////////////////////////////////
219/// Do some initializations.
220
222{
223 fRecords = nullptr;
224 fPixels = nullptr;
225 fSizes = nullptr;
226 fColumnsInfo = nullptr;
227 fNColumns = fNRows = 0;
228 fCells = nullptr;
229}
230
231////////////////////////////////////////////////////////////////////////////////
232/// Load HDU from fits file satisfying the specified filter.
233/// Returns kTRUE if success. Otherwise kFALSE.
234/// If filter == "" then the primary array is selected
235
237{
238 fitsfile *fp=nullptr;
239 int status = 0;
240 char errdescr[FLEN_STATUS+1];
241
242 // Open file with filter
243 fits_open_file(&fp, filepath_filter.Data(), READONLY, &status);
244 if (status) goto ERR;
245
246 // Read HDU number
247 int hdunum;
250
251 // Read HDU type
252 int hdutype;
253 fits_get_hdu_type(fp, &hdutype, &status);
254 if (status) goto ERR;
256
257 //Read HDU header records
258 int nkeys, morekeys;
259 char keyname[FLEN_KEYWORD+1];
260 char keyvalue[FLEN_VALUE+1];
261 char comment[FLEN_COMMENT+1];
262
263 fits_get_hdrspace(fp, &nkeys, &morekeys, &status);
264 if (status) goto ERR;
265
266 fRecords = new struct HDURecord[nkeys];
267
268 for (int i = 1; i <= nkeys; i++) {
269 fits_read_keyn(fp, i, keyname, keyvalue, comment, &status);
270 if (status) goto ERR;
271 fRecords[i-1].fKeyword = keyname;
272 fRecords[i-1].fValue = keyvalue;
273 fRecords[i-1].fComment = comment;
274 }
275
277
278 //Set extension name
279 fExtensionName = "PRIMARY"; //Default
280 for (int i = 0; i < nkeys; i++) {
281 if (fRecords[i].fKeyword == "EXTNAME") {
282 fExtensionName = fRecords[i].fValue;
283 break;
284 }
285 }
286
287 //Read HDU's data
288 if (fType == kImageHDU) {
289 Info("LoadHDU", "The selected HDU contains an Image Extension");
290
291 int param_ndims=0;
292 long *param_dimsizes;
293
294 //Read image number of dimensions
295 fits_get_img_dim(fp, &param_ndims, &status);
296 if (status) goto ERR;
297 if (param_ndims > 0) {
298 //Read image sizes in each dimension
299 param_dimsizes = new long[param_ndims];
301 if (status) {
302 delete [] param_dimsizes;
303 goto ERR;
304 }
305
308 for (int i = 0; i < param_ndims; i++) { //Use for loop to copy values instead of passing array to constructor, since 'Int_t' size may differ from 'long' size
310 }
311
312 delete [] param_dimsizes;
313
314 //Read pixels
315 int anynul;
316 long *firstpixel = new long[param_ndims];
317 double nulval = 0;
318 long npixels = 1;
319
320 for (int i = 0; i < param_ndims; i++) {
321 npixels *= (long) fSizes->GetAt(i); //Compute total number of pixels
322 firstpixel[i] = 1; //Set first pixel to read from.
323 }
324
325 double *pixels = new double[npixels];
326
328 (void *) &nulval, (void *) pixels, &anynul, &status);
329
330 if (status) {
331 delete [] firstpixel;
332 delete [] pixels;
333 goto ERR;
334 }
335
337
338 delete [] firstpixel;
339 delete [] pixels;
340
341 } else {
342 //Null array
343 fSizes = new TArrayI();
344 fPixels = new TArrayD();
345 }
346 } else {
347 Info("LoadHDU", "The selected HDU contains a Table Extension");
348
349 // Get table's number of rows and columns
350 long table_rows;
351 int table_cols;
352
353 fits_get_num_rows(fp, &table_rows, &status);
354 if (status) goto ERR;
355
357
358 fits_get_num_cols(fp, &table_cols, &status);
359 if (status) goto ERR;
360
362
363 // Allocate column info array
364 fColumnsInfo = new struct Column[table_cols];
365
366 // Read column names
367 char colname[80];
368 int colnum;
369
370 fits_get_colname(fp, CASEINSEN, (char*) "*", colname, &colnum, &status);
371 while (status == COL_NOT_UNIQUE)
372 {
373 fColumnsInfo[colnum-1].fName = colname;
374 fits_get_colname(fp, CASEINSEN, (char*) "*", colname, &colnum, &status);
375 }
376 if (status != COL_NOT_FOUND) goto ERR;
377 status = 0;
378
379 //Allocate cells
380 fCells = new union Cell [table_rows * table_cols];
381
382 // Read columns
383 int typecode;
384 long repeat, width;
386
387
388 for (colnum = 0, cellindex = 0; colnum < fNColumns; colnum++) {
389 fits_get_coltype(fp, colnum+1, &typecode, &repeat, &width, &status);
390
391 if (status) goto ERR;
392
393 if ((typecode == TDOUBLE) || (typecode == TSHORT) || (typecode == TLONG)
394 || (typecode == TFLOAT) || (typecode == TLOGICAL) || (typecode == TBIT)
395 || (typecode == TBYTE) || (typecode == TSTRING)) {
396
397 if (typecode == TSTRING) {
398 // this column contains strings
399 fColumnsInfo[colnum].fType = kString;
400
401 int dispwidth=0;
403 if (status) goto ERR;
404
405
406 char *nulval = (char*) "";
407 int anynul=0;
408 char **array;
409
410 if (dispwidth <= 0) {
411 dispwidth = 1;
412 }
413
414 array = new char* [table_rows];
415 for (long row = 0; row < table_rows; row++) {
416 array[row] = new char[dispwidth+1]; //also room for end null!
417 }
418
419 if (repeat > 0) {
420 fits_read_col(fp, TSTRING, colnum+1, 1, 1, table_rows, nulval, array, &anynul, &status);
421 if (status) {
422 for (long row = 0; row < table_rows; row++) {
423 delete [] array[row];
424 }
425 delete [] array;
426 goto ERR;
427 }
428
429 } else {
430 //No elements: set dummy
431 for (long row = 0; row < table_rows; row++) {
432 strlcpy(array[row], "-",dispwidth+1);
433 }
434 }
435
436 //Save values
437 for (long row = 0; row < table_rows; row++) {
438 fCells[cellindex++].fString = array[row];
439 }
440
441 delete [] array; //Delete temporal array holding pointer to strings, but not delete strings themselves!
442
443
444 } else {
445 // this column contains either a number or a fixed-length array
446 double nulval = 0;
447 int anynul=0;
448
450
451 double *array = nullptr;
452 char *arrayl = nullptr;
453
454 if (repeat > 0) {
455
456 if (typecode == TLOGICAL) {
457 arrayl = new char[table_rows * repeat];
459 &status);
460 if (status) {
461 delete[] arrayl;
462 goto ERR;
463 }
464 } else {
465 array = new double[table_rows * repeat]; // Hope you got a big machine! Ask China otherwise :-)
466 fits_read_col(fp, TDOUBLE, colnum + 1, 1, 1, table_rows * repeat, &nulval, array, &anynul,
467 &status);
468 if (status) {
469 delete[] array;
470 goto ERR;
471 }
472 }
473
474 } else {
475 // No elements: set dummy
476 array = new double[table_rows];
477 for (long row = 0; row < table_rows; row++) {
478 array[row] = 0.0;
479 }
480 }
481
482 // Save values
483 if (repeat == 1) {
484 // this column contains scalars
486 if (typecode == TLOGICAL) {
487 for (long row = 0; row < table_rows; row++) {
488 int temp = (signed char)arrayl[row];
489 fCells[cellindex++].fRealNumber = (double)temp;
490 }
491 delete[] arrayl;
492 } else {
493 for (long row = 0; row < table_rows; row++) {
494 fCells[cellindex++].fRealNumber = array[row];
495 }
496 delete[] array;
497 }
498 } else if (repeat > 1) {
499 // this column contains fixed-length arrays
501 if (typecode == TLOGICAL) {
502 for (long row = 0; row < table_rows; row++) {
503 double *vec = new double[repeat];
504 long offset = row * repeat;
505 for (long component = 0; component < repeat; component++) {
506 int temp = (signed char)arrayl[offset++];
507 vec[component] = (double)temp;
508 }
509 fCells[cellindex++].fRealArray = vec;
510 }
511 delete[] arrayl;
512 } else {
513 for (long row = 0; row < table_rows; row++) {
514 double *vec = new double[repeat];
515 long offset = row * repeat;
516 for (long component = 0; component < repeat; component++) {
517 vec[component] = array[offset++];
518 }
519 fCells[cellindex++].fRealArray = vec;
520 }
521 delete[] array;
522 }
523 }
524 }
525 } else if (typecode < 1) {
526 // this column contains variable-length arrays
528
529 // null variables needed by the fits_read_col
530 double nulval = 0;
531 int anynul=0;
532
533 // loop over rows to derive the total memory requirement
534 fColumnsInfo[colnum].fRowStart.assign(table_rows+1, 0);
535 fColumnsInfo[colnum].fVarLengths.assign(table_rows, 0);
536 fColumnsInfo[colnum].fRowStart[0] = 0;
537
538 for (long row = 0; row < table_rows; row++) {
539 long offset = 0;
540 long repeat_row = 0;
541
542 fits_read_descript(fp, colnum+1, row+1, &repeat_row, &offset, &status);
543 // store the starting of each row and the number of elements it contains
544 fColumnsInfo[colnum].fRowStart[row+1] = fColumnsInfo[colnum].fRowStart[row] + repeat_row;
545 fColumnsInfo[colnum].fVarLengths[row] = repeat_row;
546 }
547
548 for (long row = 0; row < table_rows; row++) {
549 // number of elements in the cell we want to read, i.e.
550 // number of elements in the variable-length array
551 int nelements = fColumnsInfo[colnum].fRowStart[row + 1] - fColumnsInfo[colnum].fRowStart[row];
552 // size of the variable-length array
553 const int size = fColumnsInfo[colnum].fVarLengths[row];
554 // vector to store the results
555 TArrayD *vec = new TArrayD(size);
556 // variable-length array have negative DATATYPE
558
559 // define the array to load the data with the CFITSIO functions
560 // a fixed arrays is needed as argument to the fits_read_col function
561 // so a vector is defined and then its `.data()` pointer is passed
562 // to fits_read_col
563 //
564 // TODO: add all type cases, for now only short, long, float and double are considered
565 //
566 if (abstype == 21) {
567 std::vector<short> data;
568 data.resize(size);
569 fits_read_col(fp, abstype, colnum + 1, row + 1, 1, nelements, &nulval, data.data(), &anynul, &status);
570 for (int i = 0; i < size; i++)
571 vec->SetAt(data[i], i);
572 } else if (abstype == 41) {
573 std::vector<int> data;
574 data.resize(size);
575 fits_read_col(fp, abstype, colnum + 1, row + 1, 1, nelements, &nulval, data.data(), &anynul, &status);
576 for (int i = 0; i < size; i++)
577 vec->SetAt(data[i], i);
578 } else if (abstype == 42) {
579 std::vector<float> data;
580 data.resize(size);
581 fits_read_col(fp, abstype, colnum + 1, row + 1, 1, nelements, &nulval, data.data(), &anynul, &status);
582 for (int i = 0; i < size; i++)
583 vec->SetAt(data[i], i);
584 } else if (abstype == 82) {
585 std::vector<double> data;
586 data.resize(size);
587 fits_read_col(fp, abstype, colnum + 1, row + 1, 1, nelements, &nulval, data.data(), &anynul, &status);
588 for (int i = 0; i < size; i++)
589 vec->SetAt(data[i], i);
590 } else {
591 Error("LoadHDU", "The variable-length array type in column %d is unknown", colnum + 1);
592 }
593 // place the vector storing the variable-length array in the corresponding cell
594 fCells[cellindex++].fRealVector = vec;
595 }
596 } else {
597 Warning("LoadHDU", "error opening FITS file. Column type %d is currently not supported", typecode);
598 }
599 }
600
601 if (hdutype == ASCII_TBL) {
602 // ASCII table
603
604 } else {
605 // Binary table
606 }
607 }
608
609 // Close file
610 fits_close_file(fp, &status);
611 return kTRUE;
612
613ERR:
615 Warning("LoadHDU", "error opening FITS file. Details: %s", errdescr);
616 status = 0;
617 if (fp) fits_close_file(fp, &status);
618 return kFALSE;
619}
620
621////////////////////////////////////////////////////////////////////////////////
622/// Get record by keyword.
623
625{
626 for (int i = 0; i < fNRecords; i++) {
627 if (fRecords[i].fKeyword == keyword) {
628 return &fRecords[i];
629 }
630 }
631 return nullptr;
632}
633
634////////////////////////////////////////////////////////////////////////////////
635/// Get the value of a given keyword. Return "" if not found.
636
638{
640 if (rec) {
641 return rec->fValue;
642 } else {
643 return *(new TString(""));
644 }
645}
646
647////////////////////////////////////////////////////////////////////////////////
648/// Print records.
649
651{
652 for (int i = 0; i < fNRecords; i++) {
653 if (fRecords[i].fComment.Length() > 0) {
654 printf("%-10s = %s / %s\n", fRecords[i].fKeyword.Data(), fRecords[i].fValue.Data(), fRecords[i].fComment.Data());
655 } else {
656 printf("%-10s = %s\n", fRecords[i].fKeyword.Data(), fRecords[i].fValue.Data());
657 }
658 }
659}
660
661////////////////////////////////////////////////////////////////////////////////
662/// Print HDU's parent file's metadata.
663
665{
666 fitsfile *fp=nullptr;
667 int status = 0;
668 char errdescr[FLEN_STATUS+1];
669 int hducount, extnum;
670 int hdutype = IMAGE_HDU;
671 const char *exttype;
672 char extname[FLEN_CARD]="PRIMARY"; //room enough
673 int verbose = (opt[0] ? 1 : 0);
674
675 // Open file with no filters: current HDU will be the primary one.
676 fits_open_file(&fp, fBaseFilePath.Data(), READONLY, &status);
677 if (status) goto ERR;
678
679 // Read HDU count
680 fits_get_num_hdus(fp, &hducount, &status);
681 if (status) goto ERR;
682 printf("Total: %d HDUs\n", hducount);
683
684 extnum = 0;
685 while(hducount) {
686 // Read HDU type
687 fits_get_hdu_type(fp, &hdutype, &status);
688 if (status) goto ERR;
689
690 if (hdutype == IMAGE_HDU) {
691 exttype="IMAGE";
692 } else if (hdutype == ASCII_TBL) {
693 exttype="ASCII TABLE";
694 } else {
695 exttype="BINARY TABLE";
696 }
697
698 //Read HDU header records
699 int nkeys, morekeys;
700 char keyname[FLEN_KEYWORD+1];
701 char keyvalue[FLEN_VALUE+1];
702 char comment[FLEN_COMMENT+1];
703
704 fits_get_hdrspace(fp, &nkeys, &morekeys, &status);
705 if (status) goto ERR;
706
707 struct HDURecord *records = new struct HDURecord[nkeys];
708
709 for (int i = 1; i <= nkeys; i++) {
710 fits_read_keyn(fp, i, keyname, keyvalue, comment, &status);
711 if (status) {
712 delete [] records;
713 goto ERR;
714 }
715
716 records[i-1].fKeyword = keyname;
717 records[i-1].fValue = keyvalue;
718 records[i-1].fComment = comment;
719
720 if (strcmp(keyname, "EXTNAME") == 0) {
721 //Extension name
723 }
724 }
725
726 //HDU info
727 printf(" [%d] %s (%s)\n", extnum, exttype, extname);
728
729 //HDU records
730 if (verbose) {
731 for (int i = 0; i < nkeys; i++) {
732 if (comment[0]) {
733 printf(" %-10s = %s / %s\n", records[i].fKeyword.Data(), records[i].fValue.Data(), records[i].fComment.Data());
734 } else {
735 printf(" %-10s = %s\n", records[i].fKeyword.Data(), records[i].fValue.Data());
736 }
737 }
738 }
739 printf("\n");
740
741 delete [] records;
742
743 hducount--;
744 extnum++;
745 if (hducount){
746 fits_movrel_hdu(fp, 1, &hdutype, &status);
747 if (status) goto ERR;
748 }
749 }
750
751 // Close file
752 fits_close_file(fp, &status);
753 return;
754
755ERR:
757 Warning("PrintFileMetadata", "error opening FITS file. Details: %s", errdescr);
758 status = 0;
759 if (fp) fits_close_file(fp, &status);
760}
761
762////////////////////////////////////////////////////////////////////////////////
763/// Print column information
764
766{
767 if (fType != kTableHDU) {
768 Warning("PrintColumnInfo", "this is not a table HDU.");
769 return;
770 }
771
772 for (Int_t i = 0; i < fNColumns; i++) {
773 switch (fColumnsInfo[i].fType){
774 case kString:
775 printf("%-20s : %s\n", fColumnsInfo[i].fName.Data(), "STRING");
776 break;
777 case kRealNumber:
778 printf("%-20s : %s\n", fColumnsInfo[i].fName.Data(), "REAL NUMBER");
779 break;
780 case kRealArray:
781 printf("%-20s : %s\n", fColumnsInfo[i].fName.Data(), "FIXED-LENGTH ARRAY");
782 break;
783 case kRealVector:
784 printf("%-20s : %s\n", fColumnsInfo[i].fName.Data(), "VARIABLE-LENGTH ARRAY");
785 break;
786 }
787 }
788}
789
790////////////////////////////////////////////////////////////////////////////////
791/// Print full table contents
792
794{
795 int printed_chars;
796
797 if (fType != kTableHDU) {
798 Warning("PrintColumnInfo", "this is not a table HDU.");
799 return;
800 }
801
802 // check that the table does not contain fixed or variable length arrays
803 // in that case is not possible to print a flattened table
804 for (Int_t col = 0; col < fNColumns; col++){
805 if (fColumnsInfo[col].fType == kRealArray) {
806 Warning("PrintColumnInfo", "The table contains column with fixed-length arrays and cannot be flattened for printing.");
807 return;
808 }
809 else if (fColumnsInfo[col].fType == kRealVector) {
810 Warning("PrintColumnInfo", "The table contains column with variable-length arrays and cannot be flattened for printing.");
811 return;
812 }
813 }
814
815 // Dump header
816 putchar('\n');
817 printed_chars = 0;
818 for (Int_t col = 0; col < fNColumns; col++) {
819 printed_chars += printf("%-10s| ", fColumnsInfo[col].fName.Data());
820 }
821 putchar('\n');
822 while(printed_chars--) {
823 putchar('-');
824 }
825 putchar('\n');
826
827 // Dump rows
828 for (Int_t row = 0; row < fNRows; row++) {
829 for (Int_t col = 0; col < fNColumns; col++) {
830 if (fColumnsInfo[col].fType == kString) {
831 printf("%-10s", fCells[col * fNRows + row].fString);
832 } else if (fColumnsInfo[col].fType == kRealNumber) {
833 printed_chars = printf("%.2lg", fCells[col * fNRows + row].fRealNumber);
834 printed_chars -= 10;
835 while (printed_chars < 0) {
836 putchar(' ');
838 }
839 }
840
841 if (col <= fNColumns - 1) printf("| ");
842 }
843 printf("\n");
844 }
845}
846
847////////////////////////////////////////////////////////////////////////////////
848/// Print metadata.
849/// Currently supported options:
850///
851/// - "" : print HDU record data
852/// - "F" : print FITS file's extension names, numbers and types
853/// - "F+": print FITS file's extension names and types and their record data
854/// - "T" : print column information when HDU is a table
855/// - "T+" : print full table (columns header and rows)
856
857void TFITSHDU::Print(const Option_t *opt) const
858{
859 if ((opt[0] == 'F') || (opt[0] == 'f')) {
860 PrintFileMetadata((opt[1] == '+') ? "+" : "");
861 } else if ((opt[0] == 'T') || (opt[0] == 't')) {
862 if (opt[1] == '+') {
863 PrintFullTable("");
864 } else {
865 PrintColumnInfo("");
866 }
867
868 } else {
870 }
871}
872
873
874////////////////////////////////////////////////////////////////////////////////
875/// Read image HDU as a displayable image. Return 0 if conversion cannot be done.
876/// If the HDU seems to be a multilayer image, 'layer' parameter can be used
877/// to retrieve the specified layer (starting from 0)
878
880{
881 if (fType != kImageHDU) {
882 Warning("ReadAsImage", "this is not an image HDU.");
883 return nullptr;
884 }
885
886 if (((fSizes->GetSize() != 2) && (fSizes->GetSize() != 3) && (fSizes->GetSize() != 4)) || ((fSizes->GetSize() == 4) && (fSizes->GetAt(3) > 1))) {
887 Warning("ReadAsImage", "could not convert image HDU to image because it has %d dimensions.", fSizes->GetSize());
888 return nullptr;
889 }
890
893
894 width = Int_t(fSizes->GetAt(0));
895 height = Int_t(fSizes->GetAt(1));
896
898
899 if ( ((fSizes->GetSize() == 2) && (layer > 0))
900 || (((fSizes->GetSize() == 3) || (fSizes->GetSize() == 4)) && (layer >= fSizes->GetAt(2)))) {
901
902 Warning("ReadAsImage", "layer out of bounds.");
903 return nullptr;
904 }
905
906 // Get the maximum and minimum pixel values in the layer to auto-stretch pixels
907 Double_t maxval = 0, minval = 0;
908 UInt_t i;
911
912 for (i = 0; i < pixels_per_layer; i++) {
914
915 if (pixvalue > maxval) {
917 }
918
919 if ((i == 0) || (pixvalue < minval)) {
921 }
922 }
923
924 //Build the image stretching pixels into a range from 0.0 to 255.0
925 //TImage *im = new TImage(width, height);
927 if (!im) return nullptr;
929
930
931 if (maxval == minval) {
932 //plain image
933 for (i = 0; i < pixels_per_layer; i++) {
934 layer_pixels->SetAt(255.0, i);
935 }
936 } else {
937 Double_t factor = 255.0 / (maxval-minval);
938 for (i = 0; i < pixels_per_layer; i++) {
940 layer_pixels->SetAt(factor * (pixvalue-minval), i) ;
941 }
942 }
943
944 if (pal == nullptr) {
945 // Default palette: grayscale palette
946 pal = new TImagePalette(256);
947 for (i = 0; i < 256; i++) {
948 pal->fPoints[i] = ((Double_t)i)/255.0;
949 pal->fColorAlpha[i] = 0xffff;
950 pal->fColorBlue[i] = pal->fColorGreen[i] = pal->fColorRed[i] = i << 8;
951 }
952 pal->fPoints[0] = 0;
953 pal->fPoints[255] = 1.0;
954
955 im->SetImage(*layer_pixels, UInt_t(width), pal);
956
957 delete pal;
958
959 } else {
960 im->SetImage(*layer_pixels, UInt_t(width), pal);
961 }
962
963 delete layer_pixels;
964
965 return im;
966}
967
968////////////////////////////////////////////////////////////////////////////////
969/// If the HDU is an image, draw the first layer of the primary array
970/// To set a title to the canvas, pass it in "opt"
971
973{
974 if (fType != kImageHDU) {
975 Warning("Draw", "cannot draw. This is not an image HDU.");
976 return;
977 }
978
979 TImage *im = ReadAsImage(0, nullptr);
980 if (im) {
984 cname.Form("%sHDU", this->GetName());
985 ctitle.Form("%d x %d", width, height);
987 im->Draw();
988 }
989}
990
991
992////////////////////////////////////////////////////////////////////////////////
993/// Read image HDU as a matrix. Return 0 if conversion cannot be done
994/// If the HDU seems to be a multilayer image, 'layer' parameter can be used
995/// to retrieve the specified layer (starting from 0) in matrix form.
996/// Options (value of 'opt'):
997/// "S": stretch pixel values to a range from 0.0 to 1.0
998
1000{
1001 if (fType != kImageHDU) {
1002 Warning("ReadAsMatrix", "this is not an image HDU.");
1003 return nullptr;
1004 }
1005
1006 if (((fSizes->GetSize() != 2) && (fSizes->GetSize() != 3) && (fSizes->GetSize() != 4)) || ((fSizes->GetSize() == 4) && (fSizes->GetAt(3) > 1))) {
1007 Warning("ReadAsMatrix", "could not convert image HDU to image because it has %d dimensions.", fSizes->GetSize());
1008 return nullptr;
1009 }
1010
1011
1012 if ( ((fSizes->GetSize() == 2) && (layer > 0))
1013 || (((fSizes->GetSize() == 3) || (fSizes->GetSize() == 4)) && (layer >= fSizes->GetAt(2)))) {
1014 Warning("ReadAsMatrix", "layer out of bounds.");
1015 return nullptr;
1016 }
1017
1020 Int_t offset;
1021 UInt_t i;
1022 TMatrixD *mat=nullptr;
1023
1024 width = Int_t(fSizes->GetAt(0));
1025 height = Int_t(fSizes->GetAt(1));
1026
1029
1030 double *layer_pixels = new double[pixels_per_layer];
1031
1032 if ((opt[0] == 'S') || (opt[0] == 's')) {
1033 //Stretch
1034 // Get the maximum and minimum pixel values in the layer to auto-stretch pixels
1035 Double_t factor, maxval=0, minval=0;
1037 for (i = 0; i < pixels_per_layer; i++) {
1038 pixvalue = fPixels->GetAt(offset + i);
1039
1040 if (pixvalue > maxval) {
1041 maxval = pixvalue;
1042 }
1043
1044 if ((i == 0) || (pixvalue < minval)) {
1045 minval = pixvalue;
1046 }
1047 }
1048
1049 if (maxval == minval) {
1050 //plain image
1051 for (i = 0; i < pixels_per_layer; i++) {
1052 layer_pixels[i] = 1.0;
1053 }
1054 } else {
1055 factor = 1.0 / (maxval-minval);
1056 mat = new TMatrixD(height, width);
1057
1058 for (i = 0; i < pixels_per_layer; i++) {
1059 layer_pixels[i] = factor * (fPixels->GetAt(offset + i) - minval);
1060 }
1061 }
1062
1063 } else {
1064 //No stretching
1065 mat = new TMatrixD(height, width);
1066
1067 for (i = 0; i < pixels_per_layer; i++) {
1068 layer_pixels[i] = fPixels->GetAt(offset + i);
1069 }
1070 }
1071
1072 if (mat) {
1073 // mat->Use(height, width, layer_pixels);
1074 memcpy(mat->GetMatrixArray(), layer_pixels, pixels_per_layer*sizeof(double));
1075 }
1076
1077 delete [] layer_pixels;
1078 return mat;
1079}
1080
1081
1082////////////////////////////////////////////////////////////////////////////////
1083/// Read image HDU as a histogram. Return 0 if conversion cannot be done.
1084/// The returned object can be TH1D, TH2D or TH3D depending on data dimensionality.
1085/// Please, check condition (returnedValue->IsA() == TH*D::%Class()) to
1086/// determine the object class.
1087///
1088/// NOTE: do not confuse with image histogram! This function interprets
1089/// the array as a histogram. It does not compute the histogram of pixel
1090/// values of an image! Here "pixels" are interpreted as number of entries.
1091
1093{
1094 if (fType != kImageHDU) {
1095 Warning("ReadAsHistogram", "this is not an image HDU.");
1096 return nullptr;
1097 }
1098
1099 TH1 *result=nullptr;
1100
1101 if ((fSizes->GetSize() != 1) && (fSizes->GetSize() != 2) && (fSizes->GetSize() != 3)) {
1102 Warning("ReadAsHistogram", "could not convert image HDU to histogram because it has %d dimensions.", fSizes->GetSize());
1103 return nullptr;
1104 }
1105
1106 if (fSizes->GetSize() == 1) {
1107 //1D
1108 UInt_t Nx = UInt_t(fSizes->GetAt(0));
1109 UInt_t x;
1110
1111 TH1D *h = new TH1D("", "", Int_t(Nx), 0, Nx-1);
1112
1113 for (x = 0; x < Nx; x++) {
1115 if (nentries < 0) nentries = 0; //Crop negative values
1116 h->Fill(x, nentries);
1117 }
1118
1119 result = h;
1120
1121 } else if (fSizes->GetSize() == 2) {
1122 //2D
1123 UInt_t Nx = UInt_t(fSizes->GetAt(0));
1124 UInt_t Ny = UInt_t(fSizes->GetAt(1));
1125 UInt_t x,y;
1126
1127 TH2D *h = new TH2D("", "", Int_t(Nx), 0, Nx-1, Int_t(Ny), 0, Ny-1);
1128
1129 for (y = 0; y < Ny; y++) {
1130 UInt_t offset = y * Nx;
1131 for (x = 0; x < Nx; x++) {
1133 if (nentries < 0) nentries = 0; //Crop negative values
1134 h->Fill(x,y, nentries);
1135 }
1136 }
1137
1138 result = h;
1139
1140 } else if (fSizes->GetSize() == 3) {
1141 //3D
1142 UInt_t Nx = UInt_t(fSizes->GetAt(0));
1143 UInt_t Ny = UInt_t(fSizes->GetAt(1));
1144 UInt_t Nz = UInt_t(fSizes->GetAt(2));
1145 UInt_t x,y,z;
1146
1147 TH3D *h = new TH3D("", "", Int_t(Nx), 0, Nx-1, Int_t(Ny), 0, Ny-1, Int_t(Nz), 0, Nz-1);
1148
1149
1150 for (z = 0; z < Nz; z++) {
1151 UInt_t offset1 = z * Nx * Ny;
1152 for (y = 0; y < Ny; y++) {
1153 UInt_t offset2 = y * Nx;
1154 for (x = 0; x < Nx; x++) {
1156 if (nentries < 0) nentries = 0; //Crop negative values
1157 h->Fill(x, y, z, nentries);
1158 }
1159 }
1160 }
1161
1162 result = h;
1163 }
1164
1165 return result;
1166}
1167
1168////////////////////////////////////////////////////////////////////////////////
1169/// Get a row from the image HDU when it's a 2D array.
1170
1172{
1173 if (fType != kImageHDU) {
1174 Warning("GetArrayRow", "this is not an image HDU.");
1175 return nullptr;
1176 }
1177
1178 if (fSizes->GetSize() != 2) {
1179 Warning("GetArrayRow", "could not get row from HDU because it has %d dimensions.", fSizes->GetSize());
1180 return nullptr;
1181 }
1182
1183 UInt_t i, offset, W,H;
1184
1185 W = UInt_t(fSizes->GetAt(0));
1186 H = UInt_t(fSizes->GetAt(1));
1187
1188
1189 if (row >= H) {
1190 Warning("GetArrayRow", "index out of bounds.");
1191 return nullptr;
1192 }
1193
1194 offset = W * row;
1195 double *v = new double[W];
1196
1197 for (i = 0; i < W; i++) {
1198 v[i] = fPixels->GetAt(offset+i);
1199 }
1200
1201 TVectorD *vec = new TVectorD(W, v);
1202
1203 delete [] v;
1204
1205 return vec;
1206}
1207
1208////////////////////////////////////////////////////////////////////////////////
1209/// Get a column from the image HDU when it's a 2D array.
1210
1212{
1213 if (fType != kImageHDU) {
1214 Warning("GetArrayColumn", "this is not an image HDU.");
1215 return nullptr;
1216 }
1217
1218 if (fSizes->GetSize() != 2) {
1219 Warning("GetArrayColumn", "could not get row from HDU because it has %d dimensions.", fSizes->GetSize());
1220 return nullptr;
1221 }
1222
1223 UInt_t i, W,H;
1224
1225 W = UInt_t(fSizes->GetAt(0));
1226 H = UInt_t(fSizes->GetAt(1));
1227
1228
1229 if (col >= W) {
1230 Warning("GetArrayColumn", "index out of bounds.");
1231 return nullptr;
1232 }
1233
1234 double *v = new double[H];
1235
1236 for (i = 0; i < H; i++) {
1237 v[i] = fPixels->GetAt(W*i+col);
1238 }
1239
1240 TVectorD *vec = new TVectorD(H, v);
1241
1242 delete [] v;
1243
1244 return vec;
1245}
1246
1247
1248////////////////////////////////////////////////////////////////////////////////
1249///Get column number given its name
1250
1252{
1253 Int_t colnum;
1254 for (colnum = 0; colnum < fNColumns; colnum++) {
1255 if (fColumnsInfo[colnum].fName == colname) {
1256 return colnum;
1257 }
1258 }
1259 return -1;
1260}
1261
1262////////////////////////////////////////////////////////////////////////////////
1263/// Get a string-typed column from a table HDU given its column index (>=0).
1264
1266{
1267 if (fType != kTableHDU) {
1268 Warning("GetTabStringColumn", "this is not a table HDU.");
1269 return nullptr;
1270 }
1271
1272 if ((colnum < 0) || (colnum >= fNColumns)) {
1273 Warning("GetTabStringColumn", "column index out of bounds.");
1274 return nullptr;
1275 }
1276
1277 if (fColumnsInfo[colnum].fType != kString) {
1278 Warning("GetTabStringColumn", "attempting to read a column that is not of type 'kString'.");
1279 return nullptr;
1280 }
1281
1283
1284 TObjArray *res = new TObjArray();
1285 for (Int_t row = 0; row < fNRows; row++) {
1286 res->Add(new TObjString(fCells[offset + row].fString));
1287 }
1288
1289 return res;
1290}
1291
1292////////////////////////////////////////////////////////////////////////////////
1293/// Get a string-typed column from a table HDU given its name
1294
1296{
1297 if (fType != kTableHDU) {
1298 Warning("GetTabStringColumn", "this is not a table HDU.");
1299 return nullptr;
1300 }
1301
1302
1304
1305 if (colnum == -1) {
1306 Warning("GetTabStringColumn", "column not found.");
1307 return nullptr;
1308 }
1309
1310 if (fColumnsInfo[colnum].fType != kString) {
1311 Warning("GetTabStringColumn", "attempting to read a column that is not of type 'kString'.");
1312 return nullptr;
1313 }
1314
1316
1317 TObjArray *res = new TObjArray();
1318 for (Int_t row = 0; row < fNRows; row++) {
1319 res->Add(new TObjString(fCells[offset + row].fString));
1320 }
1321
1322 return res;
1323}
1324
1325////////////////////////////////////////////////////////////////////////////////
1326/// Get a real number-typed column from a table HDU given its column index (>=0).
1327
1329{
1330 if (fType != kTableHDU) {
1331 Warning("GetTabRealVectorColumn", "this is not a table HDU.");
1332 return nullptr;
1333 }
1334
1335 if ((colnum < 0) || (colnum >= fNColumns)) {
1336 Warning("GetTabRealVectorColumn", "column index out of bounds.");
1337 return nullptr;
1338 }
1339
1341 Warning("GetTabRealVectorColumn", "attempting to read a column whose cells have embedded fixed-length arrays");
1342 Info("GetTabRealVectorColumn", "Use GetTabRealVectorCells() or GetTabRealVectorCell() instead.");
1343 return nullptr;
1344 } else if (fColumnsInfo[colnum].fType == kRealVector) {
1345 Warning("GetTabRealVectorColumn", "attempting to read a column whose cells have embedded variable-length arrays");
1346 Info("GetTabRealVectorColumn", "Use GetTabVarLengthCell() instead.");
1347 return nullptr;
1348 }
1349
1351
1352 Double_t *arr = new Double_t[fNRows];
1353
1354 for (Int_t row = 0; row < fNRows; row++) {
1355 arr[row] = fCells[offset + row].fRealNumber;
1356 }
1357
1358 TVectorD *res = new TVectorD();
1359 res->Use(fNRows, arr);
1360
1361 return res;
1362}
1363
1364////////////////////////////////////////////////////////////////////////////////
1365/// Get a real number-typed column from a table HDU given its name
1366
1368{
1369 if (fType != kTableHDU) {
1370 Warning("GetTabRealVectorColumn", "this is not a table HDU.");
1371 return nullptr;
1372 }
1373
1375
1376 if (colnum == -1) {
1377 Warning("GetTabRealVectorColumn", "column not found.");
1378 return nullptr;
1379 }
1380
1382 Warning("GetTabRealVectorColumn", "attempting to read a column whose cells have embedded fixed-length arrays");
1383 Info("GetTabRealVectorColumn", "Use GetTabRealVectorCells() or GetTabRealVectorCell() instead.");
1384 return nullptr;
1385 } else if (fColumnsInfo[colnum].fType == kRealVector) {
1386 Warning("GetTabRealVectorColumn", "attempting to read a column whose cells have embedded variable-length arrays");
1387 Info("GetTabRealVectorColumn", "Use GetTabVarLengthCell() instead.");
1388 return nullptr;
1389 }
1390
1392
1393 Double_t *arr = new Double_t[fNRows];
1394
1395 for (Int_t row = 0; row < fNRows; row++) {
1396 arr[row] = fCells[offset + row].fRealNumber;
1397 }
1398
1399 TVectorD *res = new TVectorD();
1400 res->Use(fNRows, arr);
1401
1402 return res;
1403}
1404
1405////////////////////////////////////////////////////////////////////////////////
1406/// Change to another HDU given by "filter".
1407/// The parameter "filter" will be appended to the
1408/// FITS file's base path. For example:
1409/// hduObject.Change("[EVENTS][TIME > 5]");
1410/// Please, see documentation of TFITSHDU(const char *filepath_with_filter) constructor
1411/// for further information.
1412
1413Bool_t TFITSHDU::Change(const char *filter)
1414{
1416 tmppath.Form("%s%s", fBaseFilePath.Data(), filter);
1417
1420
1421 if (kFALSE == LoadHDU(tmppath)) {
1422 //Failed! Restore previous hdu
1423 Warning("Change", "error changing HDU. Restoring the previous one...");
1424
1427
1428 if (kFALSE == LoadHDU(fFilePath)) {
1429 Warning("Change", "could not restore previous HDU. This object is no longer reliable!!");
1430 }
1431 return kFALSE;
1432 }
1433
1434 //Set new full path
1436 return kTRUE;
1437}
1438
1439////////////////////////////////////////////////////////////////////////////////
1440/// Change to another HDU given by extension_number
1441
1443{
1445 tmppath.Form("[%d]", extension_number);
1446
1447 return Change(tmppath.Data());
1448}
1449
1450////////////////////////////////////////////////////////////////////////////////
1451/// Get a collection of real vectors embedded in cells along a given column from a table HDU. colnum >= 0.
1452
1454{
1455 if (fType != kTableHDU) {
1456 Warning("GetTabRealVectorCells", "this is not a table HDU.");
1457 return nullptr;
1458 }
1459
1460 if ((colnum < 0) || (colnum >= fNColumns)) {
1461 Warning("GetTabRealVectorCells", "column index out of bounds.");
1462 return nullptr;
1463 }
1464
1466 Warning("GetTabRealVectorCells", "attempting to read a column whose cells have embedded variable-length arrays");
1467 Info("GetTabRealVectorCells", "Use GetTabVarLengthCell() instead.");
1468 return nullptr;
1469 }
1470
1472
1473 TObjArray *res = new TObjArray();
1474 Int_t dim = fColumnsInfo[colnum].fDim;
1475
1476 for (Int_t row = 0; row < fNRows; row++) {
1477 TVectorD *v = new TVectorD();
1478 v->Use(dim, fCells[offset + row].fRealArray);
1479 res->Add(v);
1480 }
1481
1482 //Make the collection to own the allocated TVectorD objects, so when
1483 //destroying the collection, the vectors will be destroyed too.
1484 res->SetOwner(kTRUE);
1485
1486 return res;
1487}
1488
1489////////////////////////////////////////////////////////////////////////////////
1490/// Get a collection of real vectors embedded in cells along a given column from a table HDU by name
1491
1493{
1494 if (fType != kTableHDU) {
1495 Warning("GetTabRealVectorCells", "this is not a table HDU.");
1496 return nullptr;
1497 }
1498
1500
1501 if (colnum == -1) {
1502 Warning("GetTabRealVectorCells", "column not found.");
1503 return nullptr;
1504 }
1505
1507}
1508
1509////////////////////////////////////////////////////////////////////////////////
1510/// Get a real array (with fixed or variable-length) embedded in a cell given by (row>=0, column>=0)
1511
1513{
1514 if (fType != kTableHDU) {
1515 Warning("GetTabRealVectorCell", "this is not a table HDU.");
1516 return nullptr;
1517 }
1518
1519 if ((colnum < 0) || (colnum >= fNColumns)) {
1520 Warning("GetTabRealVectorCell", "column index out of bounds.");
1521 return nullptr;
1522 }
1523
1524 if ((rownum < 0) || (rownum >= fNRows)) {
1525 Warning("GetTabRealVectorCell", "row index out of bounds.");
1526 return nullptr;
1527 }
1528
1530 Warning("GetTabRealVectorCells", "attempting to read a column whose cells have embedded variable-length arrays");
1531 Info("GetTabRealVectorCells", "Use GetTabVarLengthCell() instead.");
1532 return nullptr;
1533 }
1534
1535 TVectorD *v = new TVectorD();
1536 v->Use(fColumnsInfo[colnum].fDim, fCells[(colnum * fNRows) + rownum].fRealArray);
1537 return v;
1538}
1539
1540////////////////////////////////////////////////////////////////////////////////
1541/// Get a real vector embedded in a cell given by (row>=0, column name)
1542
1544{
1545 if (fType != kTableHDU) {
1546 Warning("GetTabRealVectorCell", "this is not a table HDU.");
1547 return nullptr;
1548 }
1549
1551
1552 if (colnum == -1) {
1553 Warning("GetTabRealVectorCell", "column not found.");
1554 return nullptr;
1555 }
1556
1558}
1559
1560////////////////////////////////////////////////////////////////////////////////
1561/// Get the name of a column given its index (column>=0).
1562/// In case of error the column name is "".
1563
1565{
1566 static TString noName;
1567 if (fType != kTableHDU) {
1568 Error("GetColumnName", "this is not a table HDU.");
1569 return noName;
1570 }
1571
1572 if ((colnum < 0) || (colnum >= fNColumns)) {
1573 Error("GetColumnName", "column index out of bounds.");
1574 return noName;
1575 }
1576 return fColumnsInfo[colnum].fName;
1577}
1578
1579////////////////////////////////////////////////////////////////////////////////
1580/// Get the variable-length array contained in a cell given by (row>=0, column name)
1581
1583
1584 if (fType != kTableHDU) {
1585 Warning("GetTabVarLengthVectorCell", "this is not a table HDU.");
1586 return nullptr;
1587 }
1588
1589 if ((colnum < 0) || (colnum >= fNColumns)) {
1590 Warning("GetTabVarLengthVectorCell", "column index out of bounds.");
1591 return nullptr;
1592 }
1593
1594 if ((rownum < 0) || (rownum >= fNRows)) {
1595 Warning("GetTabVarLengthVectorCell", "row index out of bounds.");
1596 return nullptr;
1597 }
1598
1599 return fCells[(colnum * fNRows) + rownum].fRealVector;
1600}
1601
1602////////////////////////////////////////////////////////////////////////////////
1603/// Get the variable-length array contained in a cell given by (row>=0, column name)
1604
1606{
1607 if (fType != kTableHDU) {
1608 Warning("GetTabVarLengthVectorCell", "this is not a table HDU.");
1609 return nullptr;
1610 }
1611
1613
1614 if (colnum == -1) {
1615 Warning("GetTabVarLengthVectorCell", "column not found.");
1616 return nullptr;
1617 }
1618
1620}
#define h(i)
Definition RSha256.hxx:106
size_t size(const MatrixT &matrix)
retrieve the size of a square matrix
int Int_t
Signed integer 4 bytes (int)
Definition RtypesCore.h:59
unsigned int UInt_t
Unsigned integer 4 bytes (unsigned int)
Definition RtypesCore.h:60
constexpr Bool_t kFALSE
Definition RtypesCore.h:108
double Double_t
Double 8 bytes.
Definition RtypesCore.h:73
constexpr Ssiz_t kNPOS
The equivalent of std::string::npos for the ROOT class TString.
Definition RtypesCore.h:131
constexpr Bool_t kTRUE
Definition RtypesCore.h:107
const char Option_t
Option string (const char)
Definition RtypesCore.h:80
ROOT::Detail::TRangeCast< T, true > TRangeDynCast
TRangeDynCast is an adapter class that allows the typed iteration through a TCollection.
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void data
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void char Point_t Rectangle_t WindowAttributes_t Float_t Float_t Float_t Int_t Int_t UInt_t UInt_t Rectangle_t Int_t Int_t Window_t TString Int_t GCValues_t GetPrimarySelectionOwner GetDisplay GetScreen GetColormap GetNativeEvent const char const char dpyName wid window const char font_name cursor keysym reg const char only_if_exist regb h Point_t winding char text const char depth char const char Int_t count const char ColorStruct_t color const char Pixmap_t Pixmap_t PictureAttributes_t attr const char char ret_data h unsigned char height h offset
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void char Point_t Rectangle_t WindowAttributes_t Float_t Float_t Float_t Int_t Int_t UInt_t UInt_t Rectangle_t result
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void char Point_t Rectangle_t WindowAttributes_t Float_t Float_t Float_t Int_t Int_t UInt_t UInt_t Rectangle_t Int_t Int_t Window_t TString Int_t GCValues_t GetPrimarySelectionOwner GetDisplay GetScreen GetColormap GetNativeEvent const char const char dpyName wid window const char font_name cursor keysym reg const char only_if_exist regb h Point_t winding char text const char depth char const char Int_t count const char cname
Option_t Option_t width
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void char Point_t Rectangle_t height
int nentries
TMatrixT< Double_t > TMatrixD
Definition TMatrixDfwd.h:23
TVectorT< Double_t > TVectorD
Definition TVectorDfwd.h:23
Array of doubles (64 bits per element).
Definition TArrayD.h:27
Double_t GetAt(Int_t i) const override
Definition TArrayD.h:45
Array of integers (32 bits per element).
Definition TArrayI.h:27
Double_t GetAt(Int_t i) const override
Definition TArrayI.h:45
void SetAt(Double_t v, Int_t i) override
Definition TArrayI.h:51
Int_t GetSize() const
Definition TArray.h:47
The Canvas class.
Definition TCanvas.h:23
virtual void SetOwner(Bool_t enable=kTRUE)
Set whether this collection is the owner (enable==true) of its content.
union Cell * fCells
Table cells (when fType == kTableHDU).
Definition TFITS.h:90
void PrintColumnInfo(const Option_t *) const
Print column information.
Definition TFITS.cxx:765
void _initialize_me()
Do some initializations.
Definition TFITS.cxx:221
void PrintFullTable(const Option_t *) const
Print full table contents.
Definition TFITS.cxx:793
TString & GetKeywordValue(const char *keyword)
Get the value of a given keyword. Return "" if not found.
Definition TFITS.cxx:637
@ kImageHDU
Definition TFITS.h:43
@ kTableHDU
Definition TFITS.h:44
Int_t fNRows
Number of rows (when fType == kTableHDU)
Definition TFITS.h:89
const TString & GetColumnName(Int_t colnum)
Get the name of a column given its index (column>=0).
Definition TFITS.cxx:1564
Int_t fNRecords
Number of records.
Definition TFITS.h:81
TString fFilePath
Path to HDU's file including filter.
Definition TFITS.h:78
struct Column * fColumnsInfo
Information about columns (when fType == kTableHDU)
Definition TFITS.h:87
TFITSHDU(const char *filepath_with_filter)
TFITSHDU constructor from file path with HDU selection filter.
Definition TFITS.cxx:117
Int_t fNColumns
Number of columns (when fType == kTableHDU)
Definition TFITS.h:88
TVectorD * GetTabRealVectorColumn(Int_t colnum)
Get a real number-typed column from a table HDU given its column index (>=0).
Definition TFITS.cxx:1328
~TFITSHDU() override
TFITSHDU destructor.
Definition TFITS.cxx:168
@ kString
Definition TFITS.h:48
@ kRealArray
Definition TFITS.h:50
@ kRealVector
Definition TFITS.h:51
@ kRealNumber
Definition TFITS.h:49
TArrayI * fSizes
Image sizes in each dimension (when fType == kImageHDU)
Definition TFITS.h:85
TImage * ReadAsImage(Int_t layer=0, TImagePalette *pal=nullptr)
Read image HDU as a displayable image.
Definition TFITS.cxx:879
Bool_t Change(const char *filter)
Change to another HDU given by "filter".
Definition TFITS.cxx:1413
TMatrixD * ReadAsMatrix(Int_t layer=0, Option_t *opt="")
Read image HDU as a matrix.
Definition TFITS.cxx:999
Int_t fNumber
HDU number (1=PRIMARY)
Definition TFITS.h:84
void PrintHDUMetadata(const Option_t *opt="") const
Print records.
Definition TFITS.cxx:650
struct HDURecord * fRecords
HDU metadata records.
Definition TFITS.h:80
TVectorD * GetArrayRow(UInt_t row)
Get a row from the image HDU when it's a 2D array.
Definition TFITS.cxx:1171
void Draw(Option_t *opt="") override
If the HDU is an image, draw the first layer of the primary array To set a title to the canvas,...
Definition TFITS.cxx:972
void _release_resources()
Release internal resources.
Definition TFITS.cxx:176
void PrintFileMetadata(const Option_t *opt="") const
Print HDU's parent file's metadata.
Definition TFITS.cxx:664
static void CleanFilePath(const char *filepath_with_filter, TString &dst)
Clean path from possible filter and put the result in 'dst'.
Definition TFITS.cxx:93
struct HDURecord * GetRecord(const char *keyword)
Get record by keyword.
Definition TFITS.cxx:624
void Print(const Option_t *opt="") const override
Print metadata.
Definition TFITS.cxx:857
TObjArray * GetTabStringColumn(Int_t colnum)
Get a string-typed column from a table HDU given its column index (>=0).
Definition TFITS.cxx:1265
TString fBaseFilePath
Path to HDU's file excluding filter.
Definition TFITS.h:79
TH1 * ReadAsHistogram()
Read image HDU as a histogram.
Definition TFITS.cxx:1092
TObjArray * GetTabRealVectorCells(Int_t colnum)
Get a collection of real vectors embedded in cells along a given column from a table HDU....
Definition TFITS.cxx:1453
TArrayD * fPixels
Image pixels (when fType == kImageHDU)
Definition TFITS.h:86
TVectorD * GetArrayColumn(UInt_t col)
Get a column from the image HDU when it's a 2D array.
Definition TFITS.cxx:1211
enum EHDUTypes fType
HDU type.
Definition TFITS.h:82
TVectorD * GetTabRealVectorCell(Int_t rownum, Int_t colnum)
Get a real array (with fixed or variable-length) embedded in a cell given by (row>=0,...
Definition TFITS.cxx:1512
Bool_t LoadHDU(TString &filepath_filter)
Load HDU from fits file satisfying the specified filter.
Definition TFITS.cxx:236
TString fExtensionName
Extension Name.
Definition TFITS.h:83
TArrayD * GetTabVarLengthVectorCell(Int_t rownum, Int_t colnum)
Get the variable-length array contained in a cell given by (row>=0, column name)
Definition TFITS.cxx:1582
Int_t GetColumnNumber(const char *colname)
Get column number given its name.
Definition TFITS.cxx:1251
1-D histogram with a double per channel (see TH1 documentation)
Definition TH1.h:927
TH1 is the base class of all histogram classes in ROOT.
Definition TH1.h:109
2-D histogram with a double per channel (see TH1 documentation)
Definition TH2.h:400
3-D histogram with a double per channel (see TH1 documentation)
Definition TH3.h:418
A class to define a conversion from pixel values to pixel color.
Definition TAttImage.h:33
An abstract interface to image processing library.
Definition TImage.h:29
static TImage * Create()
Create an image.
Definition TImage.cxx:34
const char * GetName() const override
Returns name of object.
Definition TNamed.h:49
TString fName
Definition TNamed.h:32
An array of TObjects.
Definition TObjArray.h:31
void Add(TObject *obj) override
Definition TObjArray.h:68
Collectable string class.
Definition TObjString.h:28
virtual void Warning(const char *method, const char *msgfmt,...) const
Issue warning message.
Definition TObject.cxx:1057
virtual void Error(const char *method, const char *msgfmt,...) const
Issue error message.
Definition TObject.cxx:1071
virtual void Info(const char *method, const char *msgfmt,...) const
Issue info message.
Definition TObject.cxx:1045
Basic string class.
Definition TString.h:138
const char * Data() const
Definition TString.h:384
@ kExact
Definition TString.h:285
void Form(const char *fmt,...)
Formats a string using a printf style format descriptor.
Definition TString.cxx:2362
TVectorT< Element > & Use(Int_t lwb, Int_t upb, Element *data)
Use the array data to fill the vector lwb..upb].
Definition TVectorT.cxx:348
Double_t y[n]
Definition legend1.C:17
Double_t x[n]
Definition legend1.C:17
#define H(x, y, z)
Short_t Abs(Short_t d)
Returns the absolute value of parameter Short_t d.
Definition TMathBase.h:124