Friday, December 18, 2009

Using the Geospatial Data Abstraction Library (GDAL) in Ubuntu

For a project, I need to process SRTM elevation data. GDAL is a software library that makes it easy to read and write raster geospatial data formats. I decided to give it a try since SRTM data uses the ESRI grid format supported by the library.

1. Install the the GDAL libraries and utilities
$ apt-get install libgdal1-1.5.0 gdal-bin
2. Sample code. The SRTM data file path is in the variable pszFilename. The code is based on the tutorial.
#include <stdio.h>
#include <gdal.h>>
int main()
{
    GDALDatasetH  hDataset;
    char *pszFilename = "/media/DATA/srtm_61_10.asc";
    GDALDriverH   hDriver;
    double        adfGeoTransform[6];
    GDALRasterBandH hBand;
    int             nBlockXSize, nBlockYSize;
    int             bGotMin, bGotMax;
    double          adfMinMax[2];
    float *pafScanline;
    int   nXSize;
    int i;

    GDALAllRegister();
    hDataset = GDALOpen( pszFilename, GA_ReadOnly );
    if( hDataset == NULL )
    {
        printf("Cannot open file.");
    }
    hDriver = GDALGetDatasetDriver( hDataset );
    printf( "Driver: %s/%s\n",
            GDALGetDriverShortName( hDriver ),
            GDALGetDriverLongName( hDriver ) );
    printf( "Size is %dx%dx%d\n",
            GDALGetRasterXSize( hDataset ), 
            GDALGetRasterYSize( hDataset ),
            GDALGetRasterCount( hDataset ) );
    if( GDALGetProjectionRef( hDataset ) != NULL )
        printf( "Projection is `%s'\n", GDALGetProjectionRef( hDataset ) );
    if( GDALGetGeoTransform( hDataset, adfGeoTransform ) == CE_None )
    {
        printf( "Origin = (%.6f,%.6f)\n",
                adfGeoTransform[0], adfGeoTransform[3] );
        printf( "Pixel Size = (%.6f,%.6f)\n",
                adfGeoTransform[1], adfGeoTransform[5] );
    }
    hBand = GDALGetRasterBand( hDataset, 1 );
    GDALGetBlockSize( hBand, &nBlockXSize, &nBlockYSize );
    printf( "Block=%dx%d Type=%s, ColorInterp=%s\n",
            nBlockXSize, nBlockYSize,
            GDALGetDataTypeName(GDALGetRasterDataType(hBand)),
            GDALGetColorInterpretationName(
            GDALGetRasterColorInterpretation(hBand)) 
    );
    adfMinMax[0] = GDALGetRasterMinimum( hBand, &bGotMin );
    adfMinMax[1] = GDALGetRasterMaximum( hBand, &bGotMax );
    if( ! (bGotMin && bGotMax) )
      GDALComputeRasterMinMax( hBand, TRUE, adfMinMax );
    printf( "Min=%.3fd, Max=%.3f\n", adfMinMax[0], adfMinMax[1] );
    if( GDALGetOverviewCount(hBand) > 0 )
      printf( "Band has %d overviews.\n", GDALGetOverviewCount(hBand));
    if( GDALGetRasterColorTable( hBand ) != NULL )
      printf( "Band has a color table with %d entries.\n", 
                 GDALGetColorEntryCount(
                 GDALGetRasterColorTable( hBand ) ) 
      );
    nXSize = GDALGetRasterBandXSize( hBand );
    pafScanline = (float *) CPLMalloc(sizeof(float)*nXSize);
    GDALRasterIO( hBand, GF_Read, 0, 0, nXSize, 1, 
                    pafScanline, nXSize, 1, GDT_Float32, 
                     0, 0 );
    for (i=0;i<20;i++){
      printf("%f\n", pafScanline[i]);
    }
}

3. Build the example
$ gcc -lgdal1.5.0 -I/usr/include/gdal -o testgdal.exe testgdal.c
4. Run
$./testgdal.exe
Driver: AAIGrid/Arc/Info ASCII Grid
Size is 6001x6001x1
Projection is `GEOGCS["WGS 84",DATUM["WGS_1984",SPHEROID["WGS 84",6378137,298.257223563,AUTHORITY["EPSG","7030"]],TOWGS84[0,0,0,0,0,0,0],AUTHORITY["EPSG","6326"]],PRIMEM["Greenwich",0,AUTHORITY["EPSG","8901"]],UNIT["degree",0.0174532925199433,AUTHORITY["EPSG","9108"]],AXIS["Lat",NORTH],AXIS["Long",EAST],AUTHORITY["EPSG","4326"]]'
Origin = (119.999583,15.000417)
Pixel Size = (0.000833,-0.000833)
Block=6001x1 Type=Int16, ColorInterp=Undefined
Min=-27.000d, Max=2375.000

Tuesday, December 15, 2009

On Research Methodologies

Last December 14, 2009 I attended the  Second CAS Student-Faculty Research Conference held at the CAS auditorium. The objective of the conference is to present the research output of the different departments and institutes of the College of Arts and Sciences. There were four plenary presentations in the conference from IMSP, DSS, IC and ICS, and IBS.  The parallel sessions held in the afternoon were divided into four clusters. 
  • Social Sciences, Humanities, MST and Education Cluster
  • Physics and Chemistry Cluster
  • Mathematics, Statistics, and Computer Science Cluster
  • Biology and Environment Cluster
I decided to attend at least one presentation from each cluster in order to know what methodologies are used by the researchers in their field and work. Based on my samples, majority of the researchers conducted some form of survey and performed statistical analysis (frequency analysis, regression, test of significance).  Conclusions were drawn from the result of the statistical analysis. Another methodology used is simulation and modeling (Kinetic-Monte-Carlo, Time Series Analysis, Neural Networks, Regression Models). The main justification for using modeling and simulation is that conducting the actual experiment is expensive and time consuming. Model validation, however, was lacking in some of the researches. Product development is another methodology employed. Researchers developed a system (hardware or software), a chemical product, or an experimental setup. The product addresses some problem in the field. A common keyword in the title of the research papers presented  is characterization. Based on my understanding, characterization is providing a quantitative(or qualitative) description of some entity, structure, substance, behaviour, problem, algorithm, etc. Thus, we can probably include characterization as one of the research methodologies. Although, characterization, I think is a broad term. Also, the validity of results of the researches are strengthened if statistical analysis was conducted.