ROOT  6.06/09
Reference Guide
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 /// 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 
65 FITS file interface class
66 
67 TFITSHDU is a class that allows extracting images and data from FITS files and contains
68 several methods to manage them.
69 */
70 
71 #include "TFITS.h"
72 #include "TROOT.h"
73 #include "TImage.h"
74 #include "TArrayI.h"
75 #include "TArrayD.h"
76 #include "TH1D.h"
77 #include "TH2D.h"
78 #include "TH3D.h"
79 #include "TVectorD.h"
80 #include "TMatrixD.h"
81 #include "TObjArray.h"
82 #include "TObjString.h"
83 #include "TCanvas.h"
84 
85 #include "fitsio.h"
86 #include <stdlib.h>
87 
89 
90 ////////////////////////////////////////////////////////////////////////////////
91 /// Clean path from possible filter and put the result in 'dst'.
92 
93 void TFITSHDU::CleanFilePath(const char *filepath_with_filter, TString &dst)
94 {
95  dst = filepath_with_filter;
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 
117 TFITSHDU::TFITSHDU(const char *filepath_with_filter)
118 {
119  _initialize_me();
120 
121  fFilePath = filepath_with_filter;
122  CleanFilePath(filepath_with_filter, fBaseFilePath);
123 
124  if (kFALSE == LoadHDU(fFilePath)) {
126  throw -1;
127  }
128 }
129 
130 ////////////////////////////////////////////////////////////////////////////////
131 /// TFITSHDU constructor from filepath and extension number.
132 
133 TFITSHDU::TFITSHDU(const char *filepath, Int_t extension_number)
134 {
135  _initialize_me();
136  CleanFilePath(filepath, fBaseFilePath);
137 
138  //Add "by extension number" filter
139  fFilePath.Form("%s[%d]", fBaseFilePath.Data(), extension_number);
140 
141  if (kFALSE == LoadHDU(fFilePath)) {
143  throw -1;
144  }
145 }
146 
147 ////////////////////////////////////////////////////////////////////////////////
148 /// TFITSHDU constructor from filepath and extension name.
149 
150 TFITSHDU::TFITSHDU(const char *filepath, const char *extension_name)
151 {
152  _initialize_me();
153  CleanFilePath(filepath, fBaseFilePath);
154 
155  //Add "by extension number" filter
156  fFilePath.Form("%s[%s]", fBaseFilePath.Data(), extension_name);
157 
158 
159  if (kFALSE == LoadHDU(fFilePath)) {
161  throw -1;
162  }
163 }
164 
165 ////////////////////////////////////////////////////////////////////////////////
166 /// TFITSHDU destructor.
167 
169 {
171 }
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 == kRealVector) {
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].fRealVector;
198  }
199  }
200  }
201 
202  delete [] fCells;
203  }
204 
205  delete [] fColumnsInfo;
206  }
207 
208 
209  }
210 }
211 
212 ////////////////////////////////////////////////////////////////////////////////
213 /// Do some initializations.
214 
216 {
217  fRecords = 0;
218  fPixels = 0;
219  fSizes = 0;
220  fColumnsInfo = 0;
221  fNColumns = fNRows = 0;
222  fCells = 0;
223 }
224 
225 ////////////////////////////////////////////////////////////////////////////////
226 /// Load HDU from fits file satisfying the specified filter.
227 /// Returns kTRUE if success. Otherwise kFALSE.
228 /// If filter == "" then the primary array is selected
229 
231 {
232  fitsfile *fp=0;
233  int status = 0;
234  char errdescr[FLEN_STATUS+1];
235 
236  // Open file with filter
237  fits_open_file(&fp, filepath_filter.Data(), READONLY, &status);
238  if (status) goto ERR;
239 
240  // Read HDU number
241  int hdunum;
242  fits_get_hdu_num(fp, &hdunum);
243  fNumber = Int_t(hdunum);
244 
245  // Read HDU type
246  int hdutype;
247  fits_get_hdu_type(fp, &hdutype, &status);
248  if (status) goto ERR;
249  fType = (hdutype == IMAGE_HDU) ? kImageHDU : kTableHDU;
250 
251  //Read HDU header records
252  int nkeys, morekeys;
253  char keyname[FLEN_KEYWORD+1];
254  char keyvalue[FLEN_VALUE+1];
255  char comment[FLEN_COMMENT+1];
256 
257  fits_get_hdrspace(fp, &nkeys, &morekeys, &status);
258  if (status) goto ERR;
259 
260  fRecords = new struct HDURecord[nkeys];
261 
262  for (int i = 1; i <= nkeys; i++) {
263  fits_read_keyn(fp, i, keyname, keyvalue, comment, &status);
264  if (status) goto ERR;
265  fRecords[i-1].fKeyword = keyname;
266  fRecords[i-1].fValue = keyvalue;
267  fRecords[i-1].fComment = comment;
268  }
269 
270  fNRecords = Int_t(nkeys);
271 
272  //Set extension name
273  fExtensionName = "PRIMARY"; //Default
274  for (int i = 0; i < nkeys; i++) {
275  if (fRecords[i].fKeyword == "EXTNAME") {
277  break;
278  }
279  }
280 
281  //Read HDU's data
282  if (fType == kImageHDU) {
283  //Image
284  int param_ndims=0;
285  long *param_dimsizes;
286 
287  //Read image number of dimensions
288  fits_get_img_dim(fp, &param_ndims, &status);
289  if (status) goto ERR;
290  if (param_ndims > 0) {
291  //Read image sizes in each dimension
292  param_dimsizes = new long[param_ndims];
293  fits_get_img_size(fp, param_ndims, param_dimsizes, &status);
294  if (status) {
295  delete [] param_dimsizes;
296  goto ERR;
297  }
298 
299  fSizes = new TArrayI(param_ndims);
300  fSizes = new TArrayI(param_ndims);
301  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
302  fSizes->SetAt(param_dimsizes[i], i);
303  }
304 
305  delete [] param_dimsizes;
306 
307  //Read pixels
308  int anynul;
309  long *firstpixel = new long[param_ndims];
310  double nulval = 0;
311  long npixels = 1;
312 
313  for (int i = 0; i < param_ndims; i++) {
314  npixels *= (long) fSizes->GetAt(i); //Compute total number of pixels
315  firstpixel[i] = 1; //Set first pixel to read from.
316  }
317 
318  double *pixels = new double[npixels];
319 
320  fits_read_pix(fp, TDOUBLE, firstpixel, npixels,
321  (void *) &nulval, (void *) pixels, &anynul, &status);
322 
323  if (status) {
324  delete [] firstpixel;
325  delete [] pixels;
326  goto ERR;
327  }
328 
329  fPixels = new TArrayD(npixels, pixels);
330 
331  delete [] firstpixel;
332  delete [] pixels;
333 
334  } else {
335  //Null array
336  fSizes = new TArrayI();
337  fPixels = new TArrayD();
338  }
339  } else {
340  // Table
341 
342  // Get table's number of rows and columns
343  long table_rows;
344  int table_cols;
345 
346  fits_get_num_rows(fp, &table_rows, &status);
347  if (status) goto ERR;
348 
349  fNRows = Int_t(table_rows);
350 
351  fits_get_num_cols(fp, &table_cols, &status);
352  if (status) goto ERR;
353 
354  fNColumns = Int_t(table_cols);
355 
356  // Allocate column info array
357  fColumnsInfo = new struct Column[table_cols];
358 
359  // Read column names
360  char colname[80];
361  int colnum;
362 
363  fits_get_colname(fp, CASEINSEN, (char*) "*", colname, &colnum, &status);
364  while (status == COL_NOT_UNIQUE)
365  {
366  fColumnsInfo[colnum-1].fName = colname;
367  fits_get_colname(fp, CASEINSEN, (char*) "*", colname, &colnum, &status);
368  }
369  if (status != COL_NOT_FOUND) goto ERR;
370  status = 0;
371 
372  //Allocate cells
373  fCells = new union Cell [table_rows * table_cols];
374 
375  // Read columns
376  int typecode;
377  long repeat, width;
378  Int_t cellindex;
379 
380 
381  for (colnum = 0, cellindex = 0; colnum < fNColumns; colnum++) {
382  fits_get_coltype(fp, colnum+1, &typecode, &repeat, &width, &status);
383  if (status) goto ERR;
384 
385  if ((typecode == TDOUBLE) || (typecode == TSHORT) || (typecode == TLONG)
386  || (typecode == TFLOAT) || (typecode == TLOGICAL) || (typecode == TBIT)
387  || (typecode == TBYTE) || (typecode == TSTRING)) {
388 
389  fColumnsInfo[colnum].fType = (typecode == TSTRING) ? kString : kRealNumber;
390 
391  if (typecode == TSTRING) {
392  // String column
393  int dispwidth=0;
394  fits_get_col_display_width(fp, colnum+1, &dispwidth, &status);
395  if (status) goto ERR;
396 
397 
398  char *nulval = (char*) "";
399  int anynul=0;
400  char **array;
401 
402  if (dispwidth <= 0) {
403  dispwidth = 1;
404  }
405 
406  array = new char* [table_rows];
407  for (long row = 0; row < table_rows; row++) {
408  array[row] = new char[dispwidth+1]; //also room for end null!
409  }
410 
411  if (repeat > 0) {
412  fits_read_col(fp, TSTRING, colnum+1, 1, 1, table_rows, nulval, array, &anynul, &status);
413  if (status) {
414  for (long row = 0; row < table_rows; row++) {
415  delete [] array[row];
416  }
417  delete [] array;
418  goto ERR;
419  }
420 
421  } else {
422  //No elements: set dummy
423  for (long row = 0; row < table_rows; row++) {
424  strlcpy(array[row], "-",dispwidth+1);
425  }
426  }
427 
428 
429  //Save values
430  for (long row = 0; row < table_rows; row++) {
431  fCells[cellindex++].fString = array[row];
432  }
433 
434  delete [] array; //Delete temporal array holding pointer to strings, but not delete strings themselves!
435 
436 
437  } else {
438  //Numeric or vector column
439  double nulval = 0;
440  int anynul=0;
441 
442  fColumnsInfo[colnum].fDim = (Int_t) repeat;
443 
444  double *array;
445 
446  if (repeat > 0) {
447  array = new double [table_rows * repeat]; //Hope you got a big machine! Ask China otherwise :-)
448  fits_read_col(fp, TDOUBLE, colnum+1, 1, 1, table_rows * repeat, &nulval, array, &anynul, &status);
449 
450  if (status) {
451  delete [] array;
452  goto ERR;
453  }
454  } else {
455  //No elements: set dummy
456  array = new double [table_rows];
457  for (long row = 0; row < table_rows; row++) {
458  array[row] = 0.0;
459  }
460  }
461 
462  //Save values
463  if (repeat == 1) {
464  //Scalar
465  for (long row = 0; row < table_rows; row++) {
466  fCells[cellindex++].fRealNumber = array[row];
467  }
468  } else if (repeat > 1) {
469  //Vector
470  for (long row = 0; row < table_rows; row++) {
471  double *vec = new double [repeat];
472  long offset = row * repeat;
473  for (long component = 0; component < repeat; component++) {
474  vec[component] = array[offset++];
475  }
476  fCells[cellindex++].fRealVector = vec;
477  }
478  }
479 
480  delete [] array;
481 
482  }
483 
484  } else {
485  Warning("LoadHDU", "error opening FITS file. Column type %d is currently not supported", typecode);
486  }
487  }
488 
489 
490 
491  if (hdutype == ASCII_TBL) {
492  //ASCII table
493 
494  } else {
495  //Binary table
496  }
497 
498  }
499 
500  // Close file
501  fits_close_file(fp, &status);
502  return kTRUE;
503 
504 ERR:
505  fits_get_errstatus(status, errdescr);
506  Warning("LoadHDU", "error opening FITS file. Details: %s", errdescr);
507  status = 0;
508  if (fp) fits_close_file(fp, &status);
509  return kFALSE;
510 }
511 
512 ////////////////////////////////////////////////////////////////////////////////
513 /// Get record by keyword.
514 
515 struct TFITSHDU::HDURecord* TFITSHDU::GetRecord(const char *keyword)
516 {
517  for (int i = 0; i < fNRecords; i++) {
518  if (fRecords[i].fKeyword == keyword) {
519  return &fRecords[i];
520  }
521  }
522  return 0;
523 }
524 
525 ////////////////////////////////////////////////////////////////////////////////
526 /// Get the value of a given keyword. Return "" if not found.
527 
528 TString& TFITSHDU::GetKeywordValue(const char *keyword)
529 {
530  HDURecord *rec = GetRecord(keyword);
531  if (rec) {
532  return rec->fValue;
533  } else {
534  return *(new TString(""));
535  }
536 }
537 
538 ////////////////////////////////////////////////////////////////////////////////
539 /// Print records.
540 
542 {
543  for (int i = 0; i < fNRecords; i++) {
544  if (fRecords[i].fComment.Length() > 0) {
545  printf("%-10s = %s / %s\n", fRecords[i].fKeyword.Data(), fRecords[i].fValue.Data(), fRecords[i].fComment.Data());
546  } else {
547  printf("%-10s = %s\n", fRecords[i].fKeyword.Data(), fRecords[i].fValue.Data());
548  }
549  }
550 }
551 
552 ////////////////////////////////////////////////////////////////////////////////
553 /// Print HDU's parent file's metadata.
554 
555 void TFITSHDU::PrintFileMetadata(const Option_t *opt) const
556 {
557  fitsfile *fp=0;
558  int status = 0;
559  char errdescr[FLEN_STATUS+1];
560  int hducount, extnum;
561  int hdutype = IMAGE_HDU;
562  const char *exttype;
563  char extname[FLEN_CARD]="PRIMARY"; //room enough
564  int verbose = (opt[0] ? 1 : 0);
565 
566  // Open file with no filters: current HDU will be the primary one.
567  fits_open_file(&fp, fBaseFilePath.Data(), READONLY, &status);
568  if (status) goto ERR;
569 
570  // Read HDU count
571  fits_get_num_hdus(fp, &hducount, &status);
572  if (status) goto ERR;
573  printf("Total: %d HDUs\n", hducount);
574 
575  extnum = 0;
576  while(hducount) {
577  // Read HDU type
578  fits_get_hdu_type(fp, &hdutype, &status);
579  if (status) goto ERR;
580 
581  if (hdutype == IMAGE_HDU) {
582  exttype="IMAGE";
583  } else if (hdutype == ASCII_TBL) {
584  exttype="ASCII TABLE";
585  } else {
586  exttype="BINARY TABLE";
587  }
588 
589  //Read HDU header records
590  int nkeys, morekeys;
591  char keyname[FLEN_KEYWORD+1];
592  char keyvalue[FLEN_VALUE+1];
593  char comment[FLEN_COMMENT+1];
594 
595  fits_get_hdrspace(fp, &nkeys, &morekeys, &status);
596  if (status) goto ERR;
597 
598  struct HDURecord *records = new struct HDURecord[nkeys];
599 
600  for (int i = 1; i <= nkeys; i++) {
601  fits_read_keyn(fp, i, keyname, keyvalue, comment, &status);
602  if (status) {
603  delete [] records;
604  goto ERR;
605  }
606 
607  records[i-1].fKeyword = keyname;
608  records[i-1].fValue = keyvalue;
609  records[i-1].fComment = comment;
610 
611  if (strcmp(keyname, "EXTNAME") == 0) {
612  //Extension name
613  strlcpy(extname, keyvalue,FLEN_CARD);
614  }
615  }
616 
617  //HDU info
618  printf(" [%d] %s (%s)\n", extnum, exttype, extname);
619 
620  //HDU records
621  if (verbose) {
622  for (int i = 0; i < nkeys; i++) {
623  if (comment[0]) {
624  printf(" %-10s = %s / %s\n", records[i].fKeyword.Data(), records[i].fValue.Data(), records[i].fComment.Data());
625  } else {
626  printf(" %-10s = %s\n", records[i].fKeyword.Data(), records[i].fValue.Data());
627  }
628  }
629  }
630  printf("\n");
631 
632  delete [] records;
633 
634  hducount--;
635  extnum++;
636  if (hducount){
637  fits_movrel_hdu(fp, 1, &hdutype, &status);
638  if (status) goto ERR;
639  }
640  }
641 
642  // Close file
643  fits_close_file(fp, &status);
644  return;
645 
646 ERR:
647  fits_get_errstatus(status, errdescr);
648  Warning("PrintFileMetadata", "error opening FITS file. Details: %s", errdescr);
649  status = 0;
650  if (fp) fits_close_file(fp, &status);
651 }
652 
653 ////////////////////////////////////////////////////////////////////////////////
654 /// Print column information
655 
657 {
658  if (fType != kTableHDU) {
659  Warning("PrintColumnInfo", "this is not a table HDU.");
660  return;
661  }
662 
663  for (Int_t i = 0; i < fNColumns; i++) {
664  printf("%-20s : %s\n", fColumnsInfo[i].fName.Data(), (fColumnsInfo[i].fType == kRealNumber) ? "REAL NUMBER" : "STRING");
665  }
666 }
667 
668 ////////////////////////////////////////////////////////////////////////////////
669 /// Print full table contents
670 
672 {
673  int printed_chars;
674 
675  if (fType != kTableHDU) {
676  Warning("PrintColumnInfo", "this is not a table HDU.");
677  return;
678  }
679 
680  // Dump header
681  putchar('\n');
682  printed_chars = 0;
683  for (Int_t col = 0; col < fNColumns; col++) {
684  printed_chars += printf("%-10s| ", fColumnsInfo[col].fName.Data());
685  }
686  putchar('\n');
687  while(printed_chars--) {
688  putchar('-');
689  }
690  putchar('\n');
691 
692  // Dump rows
693  for (Int_t row = 0; row < fNRows; row++) {
694  for (Int_t col = 0; col < fNColumns; col++) {
695  if (fColumnsInfo[col].fType == kString) {
696  printf("%-10s", fCells[col * fNRows + row].fString);
697  } else {
698  printed_chars = printf("%.2lg", fCells[col * fNRows + row].fRealNumber);
699  printed_chars -= 10;
700  while (printed_chars < 0) {
701  putchar(' ');
702  printed_chars++;
703  }
704  }
705 
706  if (col <= fNColumns - 1) printf("| ");
707  }
708  printf("\n");
709  }
710 }
711 
712 ////////////////////////////////////////////////////////////////////////////////
713 /// Print metadata.
714 /// Currently supported options:
715 ///
716 /// - "" : print HDU record data
717 /// - "F" : print FITS file's extension names, numbers and types
718 /// - "F+": print FITS file's extension names and types and their record data
719 /// - "T" : print column information when HDU is a table
720 /// - "T+" : print full table (columns header and rows)
721 
722 void TFITSHDU::Print(const Option_t *opt) const
723 {
724  if ((opt[0] == 'F') || (opt[0] == 'f')) {
725  PrintFileMetadata((opt[1] == '+') ? "+" : "");
726  } else if ((opt[0] == 'T') || (opt[0] == 't')) {
727  if (opt[1] == '+') {
728  PrintFullTable("");
729  } else {
730  PrintColumnInfo("");
731  }
732 
733  } else {
734  PrintHDUMetadata("");
735  }
736 }
737 
738 
739 ////////////////////////////////////////////////////////////////////////////////
740 /// Read image HDU as a displayable image. Return 0 if conversion cannot be done.
741 /// If the HDU seems to be a multilayer image, 'layer' parameter can be used
742 /// to retrieve the specified layer (starting from 0)
743 
745 {
746  if (fType != kImageHDU) {
747  Warning("ReadAsImage", "this is not an image HDU.");
748  return 0;
749  }
750 
751 
752  if (((fSizes->GetSize() != 2) && (fSizes->GetSize() != 3) && (fSizes->GetSize() != 4)) || ((fSizes->GetSize() == 4) && (fSizes->GetAt(3) > 1))) {
753  Warning("ReadAsImage", "could not convert image HDU to image because it has %d dimensions.", fSizes->GetSize());
754  return 0;
755  }
756 
757  Int_t width, height;
758  UInt_t pixels_per_layer;
759 
760  width = Int_t(fSizes->GetAt(0));
761  height = Int_t(fSizes->GetAt(1));
762 
763  pixels_per_layer = UInt_t(width) * UInt_t(height);
764 
765  if ( ((fSizes->GetSize() == 2) && (layer > 0))
766  || (((fSizes->GetSize() == 3) || (fSizes->GetSize() == 4)) && (layer >= fSizes->GetAt(2)))) {
767 
768  Warning("ReadAsImage", "layer out of bounds.");
769  return 0;
770  }
771 
772  // Get the maximum and minimum pixel values in the layer to auto-stretch pixels
773  Double_t maxval = 0, minval = 0;
774  UInt_t i;
775  Double_t pixvalue;
776  Int_t offset = layer * pixels_per_layer;
777 
778  for (i = 0; i < pixels_per_layer; i++) {
779  pixvalue = fPixels->GetAt(offset + i);
780 
781  if (pixvalue > maxval) {
782  maxval = pixvalue;
783  }
784 
785  if ((i == 0) || (pixvalue < minval)) {
786  minval = pixvalue;
787  }
788  }
789 
790  //Build the image stretching pixels into a range from 0.0 to 255.0
791  //TImage *im = new TImage(width, height);
792  TImage *im = TImage::Create();
793  if (!im) return 0;
794  TArrayD *layer_pixels = new TArrayD(pixels_per_layer);
795 
796 
797  if (maxval == minval) {
798  //plain image
799  for (i = 0; i < pixels_per_layer; i++) {
800  layer_pixels->SetAt(255.0, i);
801  }
802  } else {
803  Double_t factor = 255.0 / (maxval-minval);
804  for (i = 0; i < pixels_per_layer; i++) {
805  pixvalue = fPixels->GetAt(offset + i);
806  layer_pixels->SetAt(factor * (pixvalue-minval), i) ;
807  }
808  }
809 
810  if (pal == 0) {
811  // Default palette: grayscale palette
812  pal = new TImagePalette(256);
813  for (i = 0; i < 256; i++) {
814  pal->fPoints[i] = ((Double_t)i)/255.0;
815  pal->fColorAlpha[i] = 0xffff;
816  pal->fColorBlue[i] = pal->fColorGreen[i] = pal->fColorRed[i] = i << 8;
817  }
818  pal->fPoints[0] = 0;
819  pal->fPoints[255] = 1.0;
820 
821  im->SetImage(*layer_pixels, UInt_t(width), pal);
822 
823  delete pal;
824 
825  } else {
826  im->SetImage(*layer_pixels, UInt_t(width), pal);
827  }
828 
829  delete layer_pixels;
830 
831  return im;
832 }
833 
834 ////////////////////////////////////////////////////////////////////////////////
835 /// If the HDU is an image, draw the first layer of the primary array
836 /// To set a title to the canvas, pass it in "opt"
837 
839 {
840  if (fType != kImageHDU) {
841  Warning("Draw", "cannot draw. This is not an image HDU.");
842  return;
843  }
844 
845  TImage *im = ReadAsImage(0, 0);
846  if (im) {
847  Int_t width = Int_t(fSizes->GetAt(0));
848  Int_t height = Int_t(fSizes->GetAt(1));
849  TString cname, ctitle;
850  cname.Form("%sHDU", this->GetName());
851  ctitle.Form("%d x %d", width, height);
852  new TCanvas(cname, ctitle, width, height);
853  im->Draw();
854  }
855 }
856 
857 
858 ////////////////////////////////////////////////////////////////////////////////
859 /// Read image HDU as a matrix. Return 0 if conversion cannot be done
860 /// If the HDU seems to be a multilayer image, 'layer' parameter can be used
861 /// to retrieve the specified layer (starting from 0) in matrix form.
862 /// Options (value of 'opt'):
863 /// "S": stretch pixel values to a range from 0.0 to 1.0
864 
866 {
867  if (fType != kImageHDU) {
868  Warning("ReadAsMatrix", "this is not an image HDU.");
869  return 0;
870  }
871 
872  if (((fSizes->GetSize() != 2) && (fSizes->GetSize() != 3) && (fSizes->GetSize() != 4)) || ((fSizes->GetSize() == 4) && (fSizes->GetAt(3) > 1))) {
873  Warning("ReadAsMatrix", "could not convert image HDU to image because it has %d dimensions.", fSizes->GetSize());
874  return 0;
875  }
876 
877 
878  if ( ((fSizes->GetSize() == 2) && (layer > 0))
879  || (((fSizes->GetSize() == 3) || (fSizes->GetSize() == 4)) && (layer >= fSizes->GetAt(2)))) {
880  Warning("ReadAsMatrix", "layer out of bounds.");
881  return 0;
882  }
883 
884  Int_t width, height;
885  UInt_t pixels_per_layer;
886  Int_t offset;
887  UInt_t i;
888  TMatrixD *mat=0;
889 
890  width = Int_t(fSizes->GetAt(0));
891  height = Int_t(fSizes->GetAt(1));
892 
893  pixels_per_layer = UInt_t(width) * UInt_t(height);
894  offset = layer * pixels_per_layer;
895 
896  double *layer_pixels = new double[pixels_per_layer];
897 
898  if ((opt[0] == 'S') || (opt[0] == 's')) {
899  //Stretch
900  // Get the maximum and minimum pixel values in the layer to auto-stretch pixels
901  Double_t factor, maxval=0, minval=0;
902  Double_t pixvalue;
903  for (i = 0; i < pixels_per_layer; i++) {
904  pixvalue = fPixels->GetAt(offset + i);
905 
906  if (pixvalue > maxval) {
907  maxval = pixvalue;
908  }
909 
910  if ((i == 0) || (pixvalue < minval)) {
911  minval = pixvalue;
912  }
913  }
914 
915  if (maxval == minval) {
916  //plain image
917  for (i = 0; i < pixels_per_layer; i++) {
918  layer_pixels[i] = 1.0;
919  }
920  } else {
921  factor = 1.0 / (maxval-minval);
922  mat = new TMatrixD(height, width);
923 
924  for (i = 0; i < pixels_per_layer; i++) {
925  layer_pixels[i] = factor * (fPixels->GetAt(offset + i) - minval);
926  }
927  }
928 
929  } else {
930  //No stretching
931  mat = new TMatrixD(height, width);
932 
933  for (i = 0; i < pixels_per_layer; i++) {
934  layer_pixels[i] = fPixels->GetAt(offset + i);
935  }
936  }
937 
938  if (mat) {
939  // mat->Use(height, width, layer_pixels);
940  memcpy(mat->GetMatrixArray(), layer_pixels, pixels_per_layer*sizeof(double));
941  }
942 
943  delete [] layer_pixels;
944  return mat;
945 }
946 
947 
948 ////////////////////////////////////////////////////////////////////////////////
949 /// Read image HDU as a histogram. Return 0 if conversion cannot be done.
950 /// The returned object can be TH1D, TH2D or TH3D depending on data dimensionality.
951 /// Please, check condition (returnedValue->IsA() == TH*D::Class()) to
952 /// determine the object class.
953 ///
954 /// NOTE: do not confuse with image histogram! This function interprets
955 /// the array as a histogram. It does not compute the histogram of pixel
956 /// values of an image! Here "pixels" are interpreted as number of entries.
957 
959 {
960  if (fType != kImageHDU) {
961  Warning("ReadAsHistogram", "this is not an image HDU.");
962  return 0;
963  }
964 
965  TH1 *result=0;
966 
967  if ((fSizes->GetSize() != 1) && (fSizes->GetSize() != 2) && (fSizes->GetSize() != 3)) {
968  Warning("ReadAsHistogram", "could not convert image HDU to histogram because it has %d dimensions.", fSizes->GetSize());
969  return 0;
970  }
971 
972  if (fSizes->GetSize() == 1) {
973  //1D
974  UInt_t Nx = UInt_t(fSizes->GetAt(0));
975  UInt_t x;
976 
977  TH1D *h = new TH1D("", "", Int_t(Nx), 0, Nx-1);
978 
979  for (x = 0; x < Nx; x++) {
981  if (nentries < 0) nentries = 0; //Crop negative values
982  h->Fill(x, nentries);
983  }
984 
985  result = h;
986 
987  } else if (fSizes->GetSize() == 2) {
988  //2D
989  UInt_t Nx = UInt_t(fSizes->GetAt(0));
990  UInt_t Ny = UInt_t(fSizes->GetAt(1));
991  UInt_t x,y;
992 
993  TH2D *h = new TH2D("", "", Int_t(Nx), 0, Nx-1, Int_t(Ny), 0, Ny-1);
994 
995  for (y = 0; y < Ny; y++) {
996  UInt_t offset = y * Nx;
997  for (x = 0; x < Nx; x++) {
998  Int_t nentries = Int_t(fPixels->GetAt(offset + x));
999  if (nentries < 0) nentries = 0; //Crop negative values
1000  h->Fill(x,y, nentries);
1001  }
1002  }
1003 
1004  result = h;
1005 
1006  } else if (fSizes->GetSize() == 3) {
1007  //3D
1008  UInt_t Nx = UInt_t(fSizes->GetAt(0));
1009  UInt_t Ny = UInt_t(fSizes->GetAt(1));
1010  UInt_t Nz = UInt_t(fSizes->GetAt(2));
1011  UInt_t x,y,z;
1012 
1013  TH3D *h = new TH3D("", "", Int_t(Nx), 0, Nx-1, Int_t(Ny), 0, Ny-1, Int_t(Nz), 0, Nz-1);
1014 
1015 
1016  for (z = 0; z < Nz; z++) {
1017  UInt_t offset1 = z * Nx * Ny;
1018  for (y = 0; y < Ny; y++) {
1019  UInt_t offset2 = y * Nx;
1020  for (x = 0; x < Nx; x++) {
1021  Int_t nentries = Int_t(fPixels->GetAt(offset1 + offset2 + x));
1022  if (nentries < 0) nentries = 0; //Crop negative values
1023  h->Fill(x, y, z, nentries);
1024  }
1025  }
1026  }
1027 
1028  result = h;
1029  }
1030 
1031  return result;
1032 }
1033 
1034 ////////////////////////////////////////////////////////////////////////////////
1035 /// Get a row from the image HDU when it's a 2D array.
1036 
1038 {
1039  if (fType != kImageHDU) {
1040  Warning("GetArrayRow", "this is not an image HDU.");
1041  return 0;
1042  }
1043 
1044  if (fSizes->GetSize() != 2) {
1045  Warning("GetArrayRow", "could not get row from HDU because it has %d dimensions.", fSizes->GetSize());
1046  return 0;
1047  }
1048 
1049  UInt_t i, offset, W,H;
1050 
1051  W = UInt_t(fSizes->GetAt(0));
1052  H = UInt_t(fSizes->GetAt(1));
1053 
1054 
1055  if (row >= H) {
1056  Warning("GetArrayRow", "index out of bounds.");
1057  return 0;
1058  }
1059 
1060  offset = W * row;
1061  double *v = new double[W];
1062 
1063  for (i = 0; i < W; i++) {
1064  v[i] = fPixels->GetAt(offset+i);
1065  }
1066 
1067  TVectorD *vec = new TVectorD(W, v);
1068 
1069  delete [] v;
1070 
1071  return vec;
1072 }
1073 
1074 ////////////////////////////////////////////////////////////////////////////////
1075 /// Get a column from the image HDU when it's a 2D array.
1076 
1078 {
1079  if (fType != kImageHDU) {
1080  Warning("GetArrayColumn", "this is not an image HDU.");
1081  return 0;
1082  }
1083 
1084  if (fSizes->GetSize() != 2) {
1085  Warning("GetArrayColumn", "could not get row from HDU because it has %d dimensions.", fSizes->GetSize());
1086  return 0;
1087  }
1088 
1089  UInt_t i, W,H;
1090 
1091  W = UInt_t(fSizes->GetAt(0));
1092  H = UInt_t(fSizes->GetAt(1));
1093 
1094 
1095  if (col >= W) {
1096  Warning("GetArrayColumn", "index out of bounds.");
1097  return 0;
1098  }
1099 
1100  double *v = new double[H];
1101 
1102  for (i = 0; i < H; i++) {
1103  v[i] = fPixels->GetAt(W*i+col);
1104  }
1105 
1106  TVectorD *vec = new TVectorD(H, v);
1107 
1108  delete [] v;
1109 
1110  return vec;
1111 }
1112 
1113 
1114 ////////////////////////////////////////////////////////////////////////////////
1115 ///Get column number given its name
1116 
1117 Int_t TFITSHDU::GetColumnNumber(const char *colname)
1118 {
1119  Int_t colnum;
1120  for (colnum = 0; colnum < fNColumns; colnum++) {
1121  if (fColumnsInfo[colnum].fName == colname) {
1122  return colnum;
1123  }
1124  }
1125  return -1;
1126 }
1127 
1128 ////////////////////////////////////////////////////////////////////////////////
1129 /// Get a string-typed column from a table HDU given its column index (>=0).
1130 
1132 {
1133  if (fType != kTableHDU) {
1134  Warning("GetTabStringColumn", "this is not a table HDU.");
1135  return 0;
1136  }
1137 
1138  if ((colnum < 0) || (colnum >= fNColumns)) {
1139  Warning("GetTabStringColumn", "column index out of bounds.");
1140  return 0;
1141  }
1142 
1143  if (fColumnsInfo[colnum].fType != kString) {
1144  Warning("GetTabStringColumn", "attempting to read a column that is not of type 'kString'.");
1145  return 0;
1146  }
1147 
1148  Int_t offset = colnum * fNRows;
1149 
1150  TObjArray *res = new TObjArray();
1151  for (Int_t row = 0; row < fNRows; row++) {
1152  res->Add(new TObjString(fCells[offset + row].fString));
1153  }
1154 
1155  return res;
1156 }
1157 
1158 ////////////////////////////////////////////////////////////////////////////////
1159 /// Get a string-typed column from a table HDU given its name
1160 
1162 {
1163  if (fType != kTableHDU) {
1164  Warning("GetTabStringColumn", "this is not a table HDU.");
1165  return 0;
1166  }
1167 
1168 
1169  Int_t colnum = GetColumnNumber(colname);
1170 
1171  if (colnum == -1) {
1172  Warning("GetTabStringColumn", "column not found.");
1173  return 0;
1174  }
1175 
1176  if (fColumnsInfo[colnum].fType != kString) {
1177  Warning("GetTabStringColumn", "attempting to read a column that is not of type 'kString'.");
1178  return 0;
1179  }
1180 
1181  Int_t offset = colnum * fNRows;
1182 
1183  TObjArray *res = new TObjArray();
1184  for (Int_t row = 0; row < fNRows; row++) {
1185  res->Add(new TObjString(fCells[offset + row].fString));
1186  }
1187 
1188  return res;
1189 }
1190 
1191 ////////////////////////////////////////////////////////////////////////////////
1192 /// Get a real number-typed column from a table HDU given its column index (>=0).
1193 
1195 {
1196  if (fType != kTableHDU) {
1197  Warning("GetTabRealVectorColumn", "this is not a table HDU.");
1198  return 0;
1199  }
1200 
1201  if ((colnum < 0) || (colnum >= fNColumns)) {
1202  Warning("GetTabRealVectorColumn", "column index out of bounds.");
1203  return 0;
1204  }
1205 
1206  if (fColumnsInfo[colnum].fType != kRealNumber) {
1207  Warning("GetTabRealVectorColumn", "attempting to read a column that is not of type 'kRealNumber'.");
1208  return 0;
1209  } else if (fColumnsInfo[colnum].fDim > 1) {
1210  Warning("GetTabRealVectorColumn", "attempting to read a column whose cells have embedded vectors, not real scalars. Use GetTabRealVectorCells() instead.");
1211  return 0;
1212  }
1213 
1214  Int_t offset = colnum * fNRows;
1215 
1216  Double_t *arr = new Double_t[fNRows];
1217 
1218  for (Int_t row = 0; row < fNRows; row++) {
1219  arr[row] = fCells[offset + row].fRealNumber;
1220  }
1221 
1222  TVectorD *res = new TVectorD();
1223  res->Use(fNRows, arr);
1224 
1225  return res;
1226 }
1227 
1228 ////////////////////////////////////////////////////////////////////////////////
1229 /// Get a real number-typed column from a table HDU given its name
1230 
1232 {
1233  if (fType != kTableHDU) {
1234  Warning("GetTabRealVectorColumn", "this is not a table HDU.");
1235  return 0;
1236  }
1237 
1238  Int_t colnum = GetColumnNumber(colname);
1239 
1240  if (colnum == -1) {
1241  Warning("GetTabRealVectorColumn", "column not found.");
1242  return 0;
1243  }
1244 
1245  if (fColumnsInfo[colnum].fType != kRealNumber) {
1246  Warning("GetTabRealVectorColumn", "attempting to read a column that is not of type 'kRealNumber'.");
1247  return 0;
1248  } else if (fColumnsInfo[colnum].fDim > 1) {
1249  Warning("GetTabRealVectorColumn", "attempting to read a column whose cells have embedded vectors, not real scalars. Use GetTabRealVectorCells() instead.");
1250  return 0;
1251  }
1252 
1253  Int_t offset = colnum * fNRows;
1254 
1255  Double_t *arr = new Double_t[fNRows];
1256 
1257  for (Int_t row = 0; row < fNRows; row++) {
1258  arr[row] = fCells[offset + row].fRealNumber;
1259  }
1260 
1261  TVectorD *res = new TVectorD();
1262  res->Use(fNRows, arr);
1263 
1264  return res;
1265 }
1266 
1267 ////////////////////////////////////////////////////////////////////////////////
1268 /// Change to another HDU given by "filter".
1269 /// The parameter "filter" will be appended to the
1270 /// FITS file's base path. For example:
1271 /// hduObject.Change("[EVENTS][TIME > 5]");
1272 /// Please, see documentation of TFITSHDU(const char *filepath_with_filter) constructor
1273 /// for further information.
1274 
1275 Bool_t TFITSHDU::Change(const char *filter)
1276 {
1277  TString tmppath;
1278  tmppath.Form("%s%s", fBaseFilePath.Data(), filter);
1279 
1281  _initialize_me();
1282 
1283  if (kFALSE == LoadHDU(tmppath)) {
1284  //Failed! Restore previous hdu
1285  Warning("Change", "error changing HDU. Restoring the previous one...");
1286 
1288  _initialize_me();
1289 
1290  if (kFALSE == LoadHDU(fFilePath)) {
1291  Warning("Change", "could not restore previous HDU. This object is no longer reliable!!");
1292  }
1293  return kFALSE;
1294  }
1295 
1296  //Set new full path
1297  fFilePath = tmppath;
1298  return kTRUE;
1299 }
1300 
1301 ////////////////////////////////////////////////////////////////////////////////
1302 /// Change to another HDU given by extension_number
1303 
1304 Bool_t TFITSHDU::Change(Int_t extension_number)
1305 {
1306  TString tmppath;
1307  tmppath.Form("[%d]", extension_number);
1308 
1309  return Change(tmppath.Data());
1310 }
1311 
1312 ////////////////////////////////////////////////////////////////////////////////
1313 /// Get a collection of real vectors embedded in cells along a given column from a table HDU. colnum >= 0.
1314 
1316 {
1317  if (fType != kTableHDU) {
1318  Warning("GetTabRealVectorCells", "this is not a table HDU.");
1319  return 0;
1320  }
1321 
1322  if ((colnum < 0) || (colnum >= fNColumns)) {
1323  Warning("GetTabRealVectorCells", "column index out of bounds.");
1324  return 0;
1325  }
1326 
1327  if (fColumnsInfo[colnum].fType != kRealNumber) {
1328  Warning("GetTabRealVectorCells", "attempting to read a column that is not of type 'kRealNumber'.");
1329  return 0;
1330  }
1331 
1332  Int_t offset = colnum * fNRows;
1333 
1334  TObjArray *res = new TObjArray();
1335  Int_t dim = fColumnsInfo[colnum].fDim;
1336 
1337  for (Int_t row = 0; row < fNRows; row++) {
1338  TVectorD *v = new TVectorD();
1339  v->Use(dim, fCells[offset + row].fRealVector);
1340  res->Add(v);
1341  }
1342 
1343  //Make the collection to own the allocated TVectorD objects, so when
1344  //destroying the collection, the vectors will be destroyed too.
1345  res->SetOwner(kTRUE);
1346 
1347  return res;
1348 }
1349 
1350 ////////////////////////////////////////////////////////////////////////////////
1351 /// Get a collection of real vectors embedded in cells along a given column from a table HDU by name
1352 
1354 {
1355  if (fType != kTableHDU) {
1356  Warning("GetTabRealVectorCells", "this is not a table HDU.");
1357  return 0;
1358  }
1359 
1360  Int_t colnum = GetColumnNumber(colname);
1361 
1362  if (colnum == -1) {
1363  Warning("GetTabRealVectorCells", "column not found.");
1364  return 0;
1365  }
1366 
1367  return GetTabRealVectorCells(colnum);
1368 }
1369 
1370 ////////////////////////////////////////////////////////////////////////////////
1371 /// Get a real vector embedded in a cell given by (row>=0, column>=0)
1372 
1374 {
1375  if (fType != kTableHDU) {
1376  Warning("GetTabRealVectorCell", "this is not a table HDU.");
1377  return 0;
1378  }
1379 
1380  if ((colnum < 0) || (colnum >= fNColumns)) {
1381  Warning("GetTabRealVectorCell", "column index out of bounds.");
1382  return 0;
1383  }
1384 
1385  if ((rownum < 0) || (rownum >= fNRows)) {
1386  Warning("GetTabRealVectorCell", "row index out of bounds.");
1387  return 0;
1388  }
1389 
1390  if (fColumnsInfo[colnum].fType != kRealNumber) {
1391  Warning("GetTabRealVectorCell", "attempting to read a column that is not of type 'kRealNumber'.");
1392  return 0;
1393  }
1394 
1395  TVectorD *v = new TVectorD();
1396  v->Use(fColumnsInfo[colnum].fDim, fCells[(colnum * fNRows) + rownum].fRealVector);
1397 
1398  return v;
1399 }
1400 
1401 ////////////////////////////////////////////////////////////////////////////////
1402 /// Get a real vector embedded in a cell given by (row>=0, column name)
1403 
1404 TVectorD *TFITSHDU::GetTabRealVectorCell(Int_t rownum, const char *colname)
1405 {
1406  if (fType != kTableHDU) {
1407  Warning("GetTabRealVectorCell", "this is not a table HDU.");
1408  return 0;
1409  }
1410 
1411  Int_t colnum = GetColumnNumber(colname);
1412 
1413  if (colnum == -1) {
1414  Warning("GetTabRealVectorCell", "column not found.");
1415  return 0;
1416  }
1417 
1418  return GetTabRealVectorCell(rownum, colnum);
1419 }
1420 
1421 
1422 ////////////////////////////////////////////////////////////////////////////////
1423 /// Get the name of a column given its index (column>=0).
1424 /// In case of error the column name is "".
1425 
1427 {
1428  static TString noName;
1429  if (fType != kTableHDU) {
1430  Error("GetColumnName", "this is not a table HDU.");
1431  return noName;
1432  }
1433 
1434  if ((colnum < 0) || (colnum >= fNColumns)) {
1435  Error("GetColumnName", "column index out of bounds.");
1436  return noName;
1437  }
1438  return fColumnsInfo[colnum].fName;
1439 }
TString fKeyword
Definition: TFITS.h:59
Double_t * fPoints
Definition: TAttImage.h:87
virtual Int_t Fill(Double_t x)
Increment bin with abscissa X by 1.
Definition: TH1.cxx:3159
Double_t * fRealVector
Definition: TFITS.h:74
An array of TObjects.
Definition: TObjArray.h:39
TString fValue
Definition: TFITS.h:60
void SetAt(Double_t v, Int_t i)
Definition: TArrayD.h:53
THist< 2, double > TH2D
Definition: THist.h:320
Int_t fNRows
Definition: TFITS.h:89
TString fExtensionName
Definition: TFITS.h:83
Ssiz_t Length() const
Definition: TString.h:390
Collectable string class.
Definition: TObjString.h:32
TFITSHDU(const char *filepath_with_filter)
TFITSHDU constructor from file path with HDU selection filter.
Definition: TFITS.cxx:117
const char Option_t
Definition: RtypesCore.h:62
static const std::string comment("comment")
void Draw(Option_t *opt="")
If the HDU is an image, draw the first layer of the primary array To set a title to the canvas...
Definition: TFITS.cxx:838
Int_t fNRecords
Definition: TFITS.h:81
TArrayI * fSizes
Definition: TFITS.h:85
virtual void SetOwner(Bool_t enable=kTRUE)
Set whether this collection is the owner (enable==true) of its content.
TH1 * h
Definition: legend2.C:5
TMatrixD * ReadAsMatrix(Int_t layer=0, Option_t *opt="")
Read image HDU as a matrix.
Definition: TFITS.cxx:865
#define H(x, y, z)
Basic string class.
Definition: TString.h:137
TAlienJobStatus * status
Definition: TAlienJob.cxx:51
int Int_t
Definition: RtypesCore.h:41
bool Bool_t
Definition: RtypesCore.h:59
Int_t fDim
Definition: TFITS.h:67
const Bool_t kFALSE
Definition: Rtypes.h:92
virtual void Draw(Option_t *option="")
Default Draw method for all objects.
Definition: TObject.cxx:254
An abstract interface to image processing library.
Definition: TImage.h:45
void _release_resources()
Release internal resources.
Definition: TFITS.cxx:176
Bool_t LoadHDU(TString &filepath_filter)
Load HDU from fits file satisfying the specified filter.
Definition: TFITS.cxx:230
THist< 1, double > TH1D
Definition: THist.h:314
virtual void SetImage(const Double_t *, UInt_t, UInt_t, TImagePalette *=0)
Definition: TImage.h:132
Array of integers (32 bits per element).
Definition: TArrayI.h:29
UShort_t * fColorRed
Definition: TAttImage.h:88
ClassImp(TFITSHDU) void TFITSHDU
Clean path from possible filter and put the result in 'dst'.
Definition: TFITS.cxx:88
const char * Data() const
Definition: TString.h:349
~TFITSHDU()
TFITSHDU destructor.
Definition: TFITS.cxx:168
struct HDURecord * GetRecord(const char *keyword)
Get record by keyword.
Definition: TFITS.cxx:515
Double_t fRealNumber
Definition: TFITS.h:73
TVectorT< Element > & Use(Int_t lwb, Int_t upb, Element *data)
Use the array data to fill the vector lwb..upb].
Definition: TVectorT.cxx:346
Double_t x[n]
Definition: legend1.C:17
FITS file interface class.
Definition: TFITS.h:40
void _initialize_me()
Do some initializations.
Definition: TFITS.cxx:215
Double_t GetAt(Int_t i) const
Definition: TArrayD.h:47
void Print(const Option_t *opt="") const
Print metadata.
Definition: TFITS.cxx:722
TVectorD * GetArrayRow(UInt_t row)
Get a row from the image HDU when it's a 2D array.
Definition: TFITS.cxx:1037
enum EHDUTypes fType
Definition: TFITS.h:82
virtual void Error(const char *method, const char *msgfmt,...) const
Issue error message.
Definition: TObject.cxx:918
TVectorT< Double_t > TVectorD
Definition: TVectorDfwd.h:24
void PrintHDUMetadata(const Option_t *opt="") const
Print records.
Definition: TFITS.cxx:541
TObjArray * GetTabStringColumn(Int_t colnum)
Get a string-typed column from a table HDU given its column index (>=0).
Definition: TFITS.cxx:1131
TH1 * ReadAsHistogram()
Read image HDU as a histogram.
Definition: TFITS.cxx:958
TMatrixT< Double_t > TMatrixD
Definition: TMatrixDfwd.h:24
Bool_t Change(const char *filter)
Change to another HDU given by "filter".
Definition: TFITS.cxx:1275
enum EColumnTypes fType
Definition: TFITS.h:66
3-D histogram with a double per channel (see TH1 documentation)}
Definition: TH3.h:309
SVector< double, 2 > v
Definition: Dict.h:5
UShort_t * fColorAlpha
Definition: TAttImage.h:91
void Form(const char *fmt,...)
Formats a string using a printf style format descriptor.
Definition: TString.cxx:2321
unsigned int UInt_t
Definition: RtypesCore.h:42
bool verbose
virtual const char * GetName() const
Returns name of object.
Definition: TNamed.h:51
virtual const Element * GetMatrixArray() const
Definition: TMatrixT.h:211
Int_t Fill(Double_t)
Invalid Fill method.
Definition: TH3.cxx:275
UShort_t * fColorBlue
Definition: TAttImage.h:90
1-D histogram with a double per channel (see TH1 documentation)}
Definition: TH1.h:613
TString fName
Definition: TNamed.h:36
static TImage * Create()
Int_t GetSize() const
Definition: TArray.h:49
void PrintFullTable(const Option_t *) const
Print full table contents.
Definition: TFITS.cxx:671
THist< 3, double > TH3D
Definition: THist.h:326
int Ssiz_t
Definition: RtypesCore.h:63
The Canvas class.
Definition: TCanvas.h:48
TArrayD * fPixels
Definition: TFITS.h:86
void SetAt(Double_t v, Int_t i)
Definition: TArrayI.h:53
double Double_t
Definition: RtypesCore.h:55
TString fName
Definition: TFITS.h:65
TVectorD * GetTabRealVectorColumn(Int_t colnum)
Get a real number-typed column from a table HDU given its column index (>=0).
Definition: TFITS.cxx:1194
Char_t * fString
Definition: TFITS.h:72
ClassImp(TMCParticle) void TMCParticle printf(": p=(%7.3f,%7.3f,%9.3f) ;", fPx, fPy, fPz)
Int_t fNumber
Definition: TFITS.h:84
TString fComment
Definition: TFITS.h:61
int nentries
Definition: THbookFile.cxx:89
Double_t y[n]
Definition: legend1.C:17
The TH1 histogram class.
Definition: TH1.h:80
void PrintFileMetadata(const Option_t *opt="") const
Print HDU's parent file's metadata.
Definition: TFITS.cxx:555
A class to define a conversion from pixel values to pixel color.
Definition: TAttImage.h:83
Array of doubles (64 bits per element).
Definition: TArrayD.h:29
TString fBaseFilePath
Definition: TFITS.h:79
TVectorD * GetTabRealVectorCell(Int_t rownum, Int_t colnum)
Get a real vector embedded in a cell given by (row>=0, column>=0)
Definition: TFITS.cxx:1373
TString fFilePath
Definition: TFITS.h:78
typedef void((*Func_t)())
UShort_t * fColorGreen
Definition: TAttImage.h:89
static void CleanFilePath(const char *filepath_with_filter, TString &dst)
const Ssiz_t kNPOS
Definition: Rtypes.h:115
TVectorD * GetArrayColumn(UInt_t col)
Get a column from the image HDU when it's a 2D array.
Definition: TFITS.cxx:1077
union Cell * fCells
Definition: TFITS.h:90
TString & GetKeywordValue(const char *keyword)
Get the value of a given keyword. Return "" if not found.
Definition: TFITS.cxx:528
struct HDURecord * fRecords
Definition: TFITS.h:80
void Add(TObject *obj)
Definition: TObjArray.h:75
double result[121]
const TString & GetColumnName(Int_t colnum)
Get the name of a column given its index (column>=0).
Definition: TFITS.cxx:1426
Double_t GetAt(Int_t i) const
Definition: TArrayI.h:47
const Bool_t kTRUE
Definition: Rtypes.h:91
Int_t Fill(Double_t)
Invalid Fill method.
Definition: TH2.cxx:285
struct Column * fColumnsInfo
Definition: TFITS.h:87
Int_t fNColumns
Definition: TFITS.h:88
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:1315
void PrintColumnInfo(const Option_t *) const
Print column information.
Definition: TFITS.cxx:656
TImage * ReadAsImage(Int_t layer=0, TImagePalette *pal=0)
Read image HDU as a displayable image.
Definition: TFITS.cxx:744
Int_t GetColumnNumber(const char *colname)
Get column number given its name.
Definition: TFITS.cxx:1117
2-D histogram with a double per channel (see TH1 documentation)}
Definition: TH2.h:297
virtual void Warning(const char *method, const char *msgfmt,...) const
Issue warning message.
Definition: TObject.cxx:904