ROOT logo

From $ROOTSYS/tutorials/fitsio/FITS_tutorial1.C

// Open a FITS file and retrieve the first plane of the image array 
// as a TASImage object
void FITS_tutorial1()
{
   printf("\n\n--------------------------------\n");
   printf("WELCOME TO FITS tutorial #1 !!!!\n");
   printf("--------------------------------\n");
   printf("We're gonna open a FITS file that contains only the\n");
   printf("primary HDU, consisting on an image.\n");
   printf("The object you will see is a snapshot of the NGC7662 nebula,\n");
   printf("which was taken by the author on November 2009 in Barcelona (CATALONIA).\n\n");
      
   if (!gROOT->IsBatch()) {
      //printf("Press ENTER to start..."); getchar();
   }
   
   // Open primary HDU from file
   TFITSHDU *hdu = new TFITSHDU("sample1.fits");
   if (hdu == 0) {
      printf("ERROR: could not access the HDU\n"); return;
   }
   printf("File successfully open!\n");
   
   // Dump the HDUs within the FITS file
   // and also their metadata
   //printf("Press ENTER to see summary of all data stored in the file:"); getchar();
   
   hdu->Print("F+");
   
   printf("....................................\n");
   // Here we get the exposure time.
   //printf("Press ENTER to retrieve the exposure time from the HDU metadata..."); getchar();
   printf("Exposure time = %s\n", hdu->GetKeywordValue("EXPTIME").Data());

   
   // Read the primary array as a matrix,
   // selecting only layer 0.
   // This function may be useful to
   // do image processing.
   printf("....................................\n");
   printf("We can read the image as a matrix of values.\n");
   printf("This feature is useful to do image processing, e.g:\n");
   printf("histogram equalization, custom filtering, ...\n");
   //printf("Press ENTER to continue..."); getchar();
   
   TMatrixD *mat = hdu->ReadAsMatrix(0);
   mat->Print();
   delete mat;
   
   // Read the primary array as an image,
   // selecting only layer 0.
   printf("....................................\n");
   printf("Now the primary array will be read both as an image and as a histogram,\n");
   printf("and they will be shown in a canvas.\n");
   //printf("Press ENTER to continue..."); getchar();
   
   TASImage *im = hdu->ReadAsImage(0);
   
   // Read the primary array as a histogram.
   // Depending on array dimensions, returned
   // histogram will be 1D, 2D or 3D
   TH1 *hist = hdu->ReadAsHistogram();
   
   
   TCanvas *c = new TCanvas("c1", "FITS tutorial #1", 800, 300);
   c->Divide(2,1);
   c->cd(1);
   im->Draw();
   c->cd(2);
   hist->Draw("COL");

   // Clean up
   delete hdu;
}

 
 FITS_tutorial1.C:1
 FITS_tutorial1.C:2
 FITS_tutorial1.C:3
 FITS_tutorial1.C:4
 FITS_tutorial1.C:5
 FITS_tutorial1.C:6
 FITS_tutorial1.C:7
 FITS_tutorial1.C:8
 FITS_tutorial1.C:9
 FITS_tutorial1.C:10
 FITS_tutorial1.C:11
 FITS_tutorial1.C:12
 FITS_tutorial1.C:13
 FITS_tutorial1.C:14
 FITS_tutorial1.C:15
 FITS_tutorial1.C:16
 FITS_tutorial1.C:17
 FITS_tutorial1.C:18
 FITS_tutorial1.C:19
 FITS_tutorial1.C:20
 FITS_tutorial1.C:21
 FITS_tutorial1.C:22
 FITS_tutorial1.C:23
 FITS_tutorial1.C:24
 FITS_tutorial1.C:25
 FITS_tutorial1.C:26
 FITS_tutorial1.C:27
 FITS_tutorial1.C:28
 FITS_tutorial1.C:29
 FITS_tutorial1.C:30
 FITS_tutorial1.C:31
 FITS_tutorial1.C:32
 FITS_tutorial1.C:33
 FITS_tutorial1.C:34
 FITS_tutorial1.C:35
 FITS_tutorial1.C:36
 FITS_tutorial1.C:37
 FITS_tutorial1.C:38
 FITS_tutorial1.C:39
 FITS_tutorial1.C:40
 FITS_tutorial1.C:41
 FITS_tutorial1.C:42
 FITS_tutorial1.C:43
 FITS_tutorial1.C:44
 FITS_tutorial1.C:45
 FITS_tutorial1.C:46
 FITS_tutorial1.C:47
 FITS_tutorial1.C:48
 FITS_tutorial1.C:49
 FITS_tutorial1.C:50
 FITS_tutorial1.C:51
 FITS_tutorial1.C:52
 FITS_tutorial1.C:53
 FITS_tutorial1.C:54
 FITS_tutorial1.C:55
 FITS_tutorial1.C:56
 FITS_tutorial1.C:57
 FITS_tutorial1.C:58
 FITS_tutorial1.C:59
 FITS_tutorial1.C:60
 FITS_tutorial1.C:61
 FITS_tutorial1.C:62
 FITS_tutorial1.C:63
 FITS_tutorial1.C:64
 FITS_tutorial1.C:65
 FITS_tutorial1.C:66
 FITS_tutorial1.C:67
 FITS_tutorial1.C:68
 FITS_tutorial1.C:69
 FITS_tutorial1.C:70
 FITS_tutorial1.C:71
 FITS_tutorial1.C:72
 FITS_tutorial1.C:73
 FITS_tutorial1.C:74
 FITS_tutorial1.C:75
 FITS_tutorial1.C:76
 FITS_tutorial1.C:77
thumb