ROOT logo

From $ROOTSYS/tutorials/fitsio/FITS_tutorial2.C

// Open a FITS file whose primary array represents
// a spectrum (flux vs wavelength)
void FITS_tutorial2()
{
   printf("\n\n--------------------------------\n");
   printf("WELCOME TO FITS tutorial #2 !!!!\n");
   printf("--------------------------------\n");
   printf("We're gonna open a FITS file that contains the\n");
   printf("primary HDU and a little data table.\n");
   printf("The primary HDU is an array of 2 rows by 2040 columns, and\n");
   printf("they represent a radiation spectrum. The first row contains\n");
   printf("the flux data, whereas the second row the wavelengths.\n");
   printf("Data copyright: NASA\n\n");
   
   if (!gROOT->IsBatch()) {
      //printf("Press ENTER to start..."); getchar();
   }

   // Open primary HDU from file
   TFITSHDU *hdu = new TFITSHDU("sample2.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");
   printf("We are going to generate a TGraph from vectors\n");
   //printf("within the primary array. Press ENTER to continue.."); getchar();
   
   TVectorD *Y = hdu->GetArrayRow(0);
   TVectorD *X = hdu->GetArrayRow(1);
   TGraph *gr = new TGraph(*X,*Y);
      
   // Show the graphic
   TCanvas *c = new TCanvas("c1", "FITS tutorial #2", 800, 800);
   gr->Draw("BA");
   
         
   // Clean up
   delete X;
   delete Y;
   delete hdu;
}

 
 FITS_tutorial2.C:1
 FITS_tutorial2.C:2
 FITS_tutorial2.C:3
 FITS_tutorial2.C:4
 FITS_tutorial2.C:5
 FITS_tutorial2.C:6
 FITS_tutorial2.C:7
 FITS_tutorial2.C:8
 FITS_tutorial2.C:9
 FITS_tutorial2.C:10
 FITS_tutorial2.C:11
 FITS_tutorial2.C:12
 FITS_tutorial2.C:13
 FITS_tutorial2.C:14
 FITS_tutorial2.C:15
 FITS_tutorial2.C:16
 FITS_tutorial2.C:17
 FITS_tutorial2.C:18
 FITS_tutorial2.C:19
 FITS_tutorial2.C:20
 FITS_tutorial2.C:21
 FITS_tutorial2.C:22
 FITS_tutorial2.C:23
 FITS_tutorial2.C:24
 FITS_tutorial2.C:25
 FITS_tutorial2.C:26
 FITS_tutorial2.C:27
 FITS_tutorial2.C:28
 FITS_tutorial2.C:29
 FITS_tutorial2.C:30
 FITS_tutorial2.C:31
 FITS_tutorial2.C:32
 FITS_tutorial2.C:33
 FITS_tutorial2.C:34
 FITS_tutorial2.C:35
 FITS_tutorial2.C:36
 FITS_tutorial2.C:37
 FITS_tutorial2.C:38
 FITS_tutorial2.C:39
 FITS_tutorial2.C:40
 FITS_tutorial2.C:41
 FITS_tutorial2.C:42
 FITS_tutorial2.C:43
 FITS_tutorial2.C:44
 FITS_tutorial2.C:45
 FITS_tutorial2.C:46
 FITS_tutorial2.C:47
 FITS_tutorial2.C:48
 FITS_tutorial2.C:49
 FITS_tutorial2.C:50
 FITS_tutorial2.C:51
 FITS_tutorial2.C:52
thumb