GeoTIFFtoF5.c

Converts GeoTIFF files into HDF5, making use of the geotiff tag information to set the bounding box. It adds timing information to the data sets.

/* Copyright 2006 Ana Elena Buleu, CCT undergraduate researcher */

#include <stdlib.h>
#include <stdio.h>
#include <tiffio.h>
#include <string.h>
#include <hdf5.h>
#include <F5/F5F.h>
#include <F5/F5uniform.h>

/* needs libgeotiff from here http://trac.osgeo.org/geotiff/
   on FC8 run yum install libgeotiff libgeotiff-devel
 */

#include <libgeotiff/geotiff.h>
#include <libgeotiff/geo_normalize.h>
#include <libgeotiff/geovalues.h>

int GetGreyscale(TIFF *tif, char* GreyScale, uint32 stripbc);
void WriteToF5(TIFF* tif, char* GreyScale, hid_t FileID, char FileName[]);


int main(int argc, char *argv[])
{
  TIFF* tif;
  hid_t FileID;
  char FileName[50];
  int i;
  uint32 stripbc, h;
  char* GreyScale;

  FileID = H5Fcreate("CloudData.f5",H5F_ACC_TRUNC,H5P_DEFAULT,H5P_DEFAULT);    

  for (i=1; i<argc; i++)
    {
      strcpy(FileName,argv[i]);

      tif=TIFFOpen(FileName, "r");

      TIFFGetField(tif, TIFFTAG_IMAGEWIDTH, &stripbc);
     
      TIFFGetField(tif, TIFFTAG_IMAGELENGTH, &h);

      GreyScale = (char*)malloc(h*stripbc);

      GetGreyscale(tif, GreyScale, stripbc);

      WriteToF5(tif, GreyScale, FileID, FileName);

      free(GreyScale);

      TIFFClose(tif);
    }

  F5Xclose(FileID);

  return 0;
}

int GetGreyscale(TIFF *tif, char* GreyScale, uint32 stripbc)
{
  tstrip_t strip;
  tsize_t scanline;
  tdata_t buf;
  uint32 row, h, nrow;
  uint32 rowsperstrip = (uint32)-1;
  tsample_t s, samplesperpixel;
  int i, j;


  tstrip_t nstrips = TIFFNumberOfStrips(tif);
  const char* what = TIFFIsTiled(tif) ? "Tile" : "Strip";

  printf("Tiff is Tiled: %s\n\n", what);

  scanline = TIFFScanlineSize(tif);

  buf = (unsigned char *)_TIFFmalloc(TIFFScanlineSize(tif));
  if (buf) 
    {
      TIFFGetField(tif, TIFFTAG_IMAGELENGTH, &h);
      TIFFGetField(tif, TIFFTAG_ROWSPERSTRIP, &rowsperstrip);
      TIFFGetField(tif, TIFFTAG_SAMPLESPERPIXEL, &samplesperpixel);

      for (row = 0; row < h; row += rowsperstrip) 
        {
          for (s = 0; s < samplesperpixel; s++) 
            {
              nrow = (row+rowsperstrip > h ? h-row : rowsperstrip);
                strip = TIFFComputeStrip(tif, row, s);

              if (TIFFReadEncodedStrip(tif, strip, buf, nrow*scanline) < 0) 
          
                return 1;
          else bcopy((unsigned char*)buf, &GreyScale[row*stripbc], stripbc);  
              }
        }
      _TIFFfree(buf);
    }
  else return 1;

  return 0;
}
     
void WriteToF5(TIFF* tif, char* GreyScale, hid_t FileID, char FileName[])
{
  int i, j;
  F5Path* fpath;
  hsize_t dims[3]; /* Number of Values in each dimension */
  F5_vec3_point_t Origin; /* Starting Point of the Volume */
  F5_vec3_float_t Delta; /* Stepsize between 2 Points on the Grid */
  double Time = 0; /* Used to store the time value */
  uint32 w, l;
  char D[3], H[3], M[3];
  GTIF* gtif;
  double ulx, uly, llx, lly, lrx, lry;


  D[0]=FileName[7];
  D[1]=FileName[8];
  D[2]='\0';
  H[0]=FileName[10];
  H[1]=FileName[11];
  H[2]='\0';
  M[0]=FileName[12];
  M[1]=FileName[13];
  M[2]='\0';
  Time=60*24*(atoi(D)-27)+atoi(H)*60+atoi(M);

  TIFFGetField(tif, TIFFTAG_IMAGELENGTH, &l);
  TIFFGetField(tif,TIFFTAG_IMAGEWIDTH, &w);

  ulx=0.0;
  uly=0.0;
  llx=0.0;
  lly=(double)l;
  lrx=(double)w;
  lry=(double)l;

  gtif=GTIFNew(tif);

  if(!GTIFImageToPCS(gtif, &ulx, &uly)) return;
  if(!GTIFImageToPCS(gtif, &llx, &lly)) return;
  if(!GTIFImageToPCS(gtif, &lrx, &lry)) return;

  GTIFFree(gtif);

  dims[0] = w; /* Dimensions */
  dims[1] = l;
  dims[2] = 1;
  Origin.x = ulx;              /* Volume starts here */
  Origin.y = uly;
  Origin.z = 0; 
  Delta.x = (lrx-llx)/(int)dims[0];      /* Stepsize on Grid */
  Delta.y = (lly-uly)/(int)dims[1];
  Delta.z = 2.0/(int)dims[2];

                /* Write the Data */
  fpath = F5Fwrite_uniform_cartesian3D(
                                        FileID,           /* File identifier given by HDF5 */
                                        Time,             /* Time Value of the field */
                                        "Clouds",         /* Name of the Grid */
                                        &Origin,          /* Starting Point of the Volume */
                                        &Delta,           /* Spacing between Grid points */
                                        dims,             /* Number of Values in each direction */                                        "Infrared",      /* Name of the Field */
                                        H5T_NATIVE_UCHAR, /* Type of the Values in the Field */
                                        GreyScale,       /* Pointer to the Data */
                                        NULL,             /* Coordinate System - standard chart*/
                                        F5P_DEFAULT);     /* Default dataset properties */
  F5close(fpath); /* Close the DataSet */
}