root/MGET/Branches/Jason/VisualStudioSolutions/vectorshp/vectorshp.cpp @ 905

// vectorshp.cpp : Defines the entry point for the console application.
//

#include <iostream>
#include <tchar.h>
#include <math.h>
#include <windows.h>

#include "shapefil.h"

const DWORD BUF_SIZE = 256*1024;

enum DataType { Int8_DataType, Uint8_DataType, Int16_DataType, Uint16_DataType, Int32_DataType, Float_DataType, Double_DataType, };

const int DataTypeSizes[] = {1, 1, 2, 2, 4, 4, 8 };

static void PrintUsage()
{
    // Do not extend strings past column 91, so they come out properly
    // in an 80-column CMD shell window.

    printf("NAME\n");
    printf("\n");
    printf("     vectorshp - creates an ESRI shapefile of vectors from two 32-bit floating\n");
    printf("                 point binary rasters that represent x and y components\n");
    printf("\n");
    printf("SYNOPSIS\n");
    printf("\n");
    printf("     vectorshp xbinary ybinary outfile datatype cols rows cellsize xllcorner\n");
    printf("               yllcorner [options...]\n");
    printf("\n");
    printf("     xbinary     Binary raster representing the x components\n");
    printf("     ybinary     Binary raster representing the y components\n");
    printf("     outfile     Shapefile to create (must not exist already)\n");
    printf("     datatype    Data type of the binary rasters, one of: int8, uint8, int16,\n");
    printf("                 uint16, int32, float, or double\n");
    printf("     cols        Number of columns in the binary files\n");
    printf("     rows        Number of rows in the binary files\n");
    printf("     cellsize    Cell size of the binary files\n");
    printf("     xllcorner   X coordinate of the lower left corner of the binary files\n");
    printf("     yllcorner   Y coordinate of the lower left corner of the binary files\n");
    printf("\n");
    printf("DESCRIPTION\n");
    printf("\n");
    printf("     vectorshp accepts as input two binary rasters having the same data type\n");
    printf("     and dimensions. The first raster represents the X component of vectors\n");
    printf("     and the second represents the Y component. vectorshp outputs an ESRI\n");
    printf("     shapefile containing lines representing the vectors. The shapefile is\n");
    printf("     useful for visualizing vector fields, like a quiver plot in matlab. In\n");
    printf("     oceanography, vectorshp has been used to visualize geostrophic currents\n");
    printf("     estimated from sea surface height data and sea surface winds estimated by\n");
    printf("     the QuickSCAT satellite.\n");
    printf("\n");
    printf("     The output shapefile contains two attributes, Magnitude and Direction,\n");
    printf("     which describe the vector in a polar coordinate system. Magnitude is in\n");
    printf("     the units of the original X and Y components. Direction is in degrees,\n");
    printf("     with 0 being due east, +90 due north, -90 due south, and +180 being due\n");
    printf("     west.\n");
    printf("\n");
    printf("     The lines originate at the coordinates of the centers of the raster\n");
    printf("     cells, allowing you to overlay the shapefile on the original rasters.\n");
    printf("\n");
    printf("     A special situation arises when your rasters are in a geographic\n");
    printf("     projection but the vector components are in something other than degrees.\n");
    printf("     For example, Aviso produces geostrophic current rasters projected in the\n");
    printf("     Topex/Poseidon geographic projection with 0.25 degree cell size. The\n");
    printf("     cells give current velocity in cm/s. You must consider three things when\n");
    printf("     you create vector shapefiles.\n");
    printf("\n");
    printf("     First, the length of shapefile lines are computed as cm/s using simple\n");
    printf("     Pythagorean math. They are not computed as degrees/s, even though the\n");
    printf("     shapefile is in a geographic projection, not a cartesian projection\n");
    printf("     such as Albers Equal Area. As a result, a 2 cm/s vector will span 2\n");
    printf("     degrees. The current is clearly not moving at 2 degrees/s. If you want\n");
    printf("     the lines to be scaled according to the projection, you must first\n");
    printf("     project your input rasters to a cartesian system that uses the same\n");
    printf("     units as the vector components (e.g. Albers Equal Area with cm units).\n");
    printf("\n");
    printf("     Second, depending on the cell size of and units of your rasters, the\n");
    printf("     shapefile lines may be excessively long, resulting in a scrambled mess.\n");
    printf("     Use the -s option to reduce the length of the lines. For example, for\n");
    printf("     Aviso geostrophic currents, I typically use -s 0.01.\n");
    printf("\n");
    printf("     Finally, vectorshp does not define the ArcGIS projection of the\n");
    printf("     shapefile. You must do this yourself (or wrap vectorshp in a\n");
    printf("     geoprocessing script, as is done by Marine Geospatial Ecology Tools).\n");
    printf("     After defining the projection and loading the shapefiles into ArcMap, you\n");
    printf("     will get a warning saying the shapefile extends beyond the normal extent\n");
    printf("     for geographic projections (+/-180 longitude, +/- 90 latitude). This is\n");
    printf("     not a problem; it is expected.\n");
    printf("\n");
    printf("OPTIONS\n");
    printf("\n");
    printf("     -n val      Specifies that val is the \"NODATA\" value. Lines will not be\n");
    printf("                 created if either the X or Y component has this value.\n");
    printf("\n");
    printf("     -s fac      Specifies that lines should be scaled by fac (for example,\n");
    printf("                 -s 10 will create lines 10 times longer than normal). Use\n");
    printf("                 this option to create vectors that can be overlaid on maps\n");
    printf("                 without appearing too cluttered or too small to be visible.\n");
    printf("                 The -s option does not affect the Magnitude attribute of the\n");
    printf("                 lines; that attribute will still contain the unscaled values.\n");
    printf("\n");
    printf("     -u          Specifies that lengths of the lines should not be\n");
    printf("                 proportional to the magnitude; instead, all lines will have\n");
    printf("                 the same length. By default, this will be the cell size. If\n");
    printf("                 the -s parameter is also specified, it scales this length\n");
    printf("                 (e.g. -s 0.5 will create lines with a length of 1/2 of the\n");
    printf("                 cell size. The -u parameter does not affect the Magnitude\n");
    printf("                 attribute of the lines; that attribute will still contain the.\n");
    printf("                 unscaled magnitude of each line.\n");
    printf("\n");
    printf("RELEASE HISTORY\n");
    printf("\n");
    printf("     1.1.0 (3-Feb-2012) - Added -u option.\n");
    printf("     1.0.0 (2-Aug-2006) - Initial version.\n");
    printf("\n");
    printf("COPYRIGHT AND LICENSE\n");
    printf("\n");
    printf("     Copyright (c) 2012 Jason J. Roberts.\n");
    printf("\n");
    printf("     Compiled with portions of shapelib 1.2.10 Copyright (c) 1999, Frank\n");
    printf("     Warmerdam.\n");
    printf("\n");
    printf("     Thanks Frank, and to your contributors, for a great library!\n");
    printf("\n");
    printf("     Permission is hereby granted, free of charge, to any person obtaining a\n");
    printf("     copy of this software (the \"Software\"), to deal in the Software without\n");
    printf("     restriction, including without limitation the rights to use, copy, modify,\n");
    printf("     merge, publish, distribute, sublicense, and/or sell copies of the\n");
    printf("     Software, and to permit persons to whom the Software is furnished to do\n");
    printf("     so, subject to the following conditions:\n");
    printf("\n");
    printf("     The above copyright notice and this permission notice shall be included in\n");
    printf("     all copies or substantial portions of the Software.\n");
    printf("\n");
    printf("     THE SOFTWARE IS PROVIDED \"AS IS\", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR\n");
    printf("     IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,\n");
    printf("     FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL\n");
    printf("     THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER\n");
    printf("     LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING\n");
    printf("     FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER\n");
    printf("     DEALINGS IN THE SOFTWARE.\n");
    printf("\n");
    printf("CONTACT INFORMATION\n");
    printf("\n");
    printf("     Please email questions or feedback to jason.roberts@duke.edu.\n");
}

static bool ParseInteger(const char *pszStr, int *pValue)
{
    if (*pszStr == '\0')
        return 0;

    char *pEnd = "z";
    int i = strtol(pszStr, &pEnd, 0);

    if (*pEnd != '\0' || (i == LONG_MAX || i == LONG_MIN) && errno == ERANGE)
        return false;

    *pValue = i;

    return true;
}

static bool ParseDouble(const char *pszStr, double *pValue)
{
    if (*pszStr == '\0')
        return 0;

    // Parse the value.

    char *pEnd = "z";
    double f = strtod(pszStr, &pEnd);

    if (*pEnd != '\0' || (f == HUGE_VAL || f == -HUGE_VAL) && errno == ERANGE)
        return false;

    *pValue = f;

    return true;
}

static bool ParseFloat(const char *pszStr, float *pValue)
{
    double d;

    bool retval = ParseDouble(pszStr, &d);

    *pValue = (float) d;

    return retval;
}

static char *ParseCommandLineArguments(int argc, char *argv[], char **pXBinary, char **pYBinary, char **pOutFile, DataType *pDataType, int *pCols, int *pRows, double *pCellSize, double *pXLowerLeftCorner, double *pYLowerLeftCorner, int *pUseNoData, int *piNoDataValue, double *pdNoDataValue, double *pScaleFactor, int *pUniformLines)
{
    // Perform basic range checks.

    if (argc < 10)
        return "ERROR: Too few command line arguments specified.\n";

    if (argc > 15)
        return "ERROR: Too many command line arguments specified.\n";

    // If the files are invalid, we will detect it when we fail to open them.

    *pXBinary = argv[1];
    *pYBinary = argv[2];
    *pOutFile = argv[3];

    // Extract the remaining arguments.

    int iArgIndex = 4;

    if (_stricmp(argv[iArgIndex], "int8") == 0)
        *pDataType = Int8_DataType;
    else if (_stricmp(argv[iArgIndex], "uint8") == 0)
        *pDataType = Uint8_DataType;
    else if (_stricmp(argv[iArgIndex], "int16") == 0)
        *pDataType = Int16_DataType;
    else if (_stricmp(argv[iArgIndex], "uint16") == 0)
        *pDataType = Uint16_DataType;
    else if (_stricmp(argv[iArgIndex], "int32") == 0)
        *pDataType = Int32_DataType;
    else if (_stricmp(argv[iArgIndex], "float") == 0)
        *pDataType = Float_DataType;
    else if (_stricmp(argv[iArgIndex], "double") == 0)
        *pDataType = Double_DataType;
    else
        return "ERROR: An invalid data type was specified. The valid data types are int8, uint8, int16, uint16, int32, float and double.";
    iArgIndex++;

    if (!ParseInteger(argv[iArgIndex], pCols) || *pCols < 1)
        return "ERROR: An invalid number of columns was specified. The number of columns must be an integer greater than or equal to 1 and less than or equal to 2147483647.";
    iArgIndex++;

    if (!ParseInteger(argv[iArgIndex], pRows) || *pRows < 1)
        return "ERROR: An invalid number of rows was specified. The number of rows must be an integer greater than or equal to 1 and less than or equal to 2147483647.";
    iArgIndex++;

    if (!ParseDouble(argv[iArgIndex], pCellSize) || *pCellSize < 0.0)
        return "ERROR: An invalid cell size was specified. The cell size must be the decimal representation of a double-precision (64-bit) floating point number greater than zero.";
    iArgIndex++;

    if (!ParseDouble(argv[iArgIndex], pXLowerLeftCorner))
        return "ERROR: An invalid X coordinate of the lower-left corner was specified. The coordinate must be the decimal representation of a double-precision (64-bit) floating point number.";
    iArgIndex++;

    if (!ParseDouble(argv[iArgIndex], pYLowerLeftCorner))
        return "ERROR: An invalid X coordinate of the lower-left corner was specified. The coordinate must be the decimal representation of a double-precision (64-bit) floating point number.";
    iArgIndex++;

    *pUseNoData = 0;
    *pScaleFactor = 1.0;

    while (iArgIndex < argc)
    {
        if (argv[iArgIndex][0] != '-' || argv[iArgIndex][1] == 0 || argv[iArgIndex][2] != 0 )
            return "ERROR: Unrecognized command line argument.\n";

        switch (argv[iArgIndex][1])
        {
            case 'n':
            case 'N':
                iArgIndex++;
                if (*pDataType == Int8_DataType || *pDataType == Uint8_DataType || *pDataType == Int16_DataType || *pDataType == Uint16_DataType || *pDataType == Int32_DataType)
                {
                    if (!ParseInteger(argv[iArgIndex], piNoDataValue))
                        return "ERROR: An invalid NODATA value was specified. The NODATA value must be an integer.";
                    if (*pDataType == Int8_DataType && (*piNoDataValue < -128 || *piNoDataValue > 127))
                        return "ERROR: An invalid NODATA value was specified. The NODATA value must be an integer between -128 and 127, inclusive.";
                    if (*pDataType == Uint8_DataType && (*piNoDataValue < 0 || *piNoDataValue > 255))
                        return "ERROR: An invalid NODATA value was specified. The NODATA value must be an integer between 0 and 255, inclusive.";
                    if (*pDataType == Int16_DataType && (*piNoDataValue < -32768 || *piNoDataValue > 32767))
                        return "ERROR: An invalid NODATA value was specified. The NODATA value must be an integer between -32768 and 32767, inclusive.";
                    if (*pDataType == Uint16_DataType && (*piNoDataValue < 0 || *piNoDataValue > 65535))
                        return "ERROR: An invalid NODATA value was specified. The NODATA value must be an integer between 0 and 65535, inclusive.";
                    if (*pDataType == Int32_DataType && (*piNoDataValue < (0-2147483648) || *piNoDataValue > 2147483647))
                        return "ERROR: An invalid NODATA value was specified. The NODATA value must be an integer between -2147483648 and 2147483647, inclusive.";
                }
                else
                {
                    if (!ParseDouble(argv[iArgIndex], pdNoDataValue))
                        return "ERROR: An invalid NODATA value was specified. The NODATA value must be a decimal number.";
                }
                *pUseNoData = 1;
                iArgIndex++;
                break;

            case 's':
            case 'S':
                iArgIndex++;
                if (!ParseDouble(argv[iArgIndex], pScaleFactor) || *pScaleFactor <= 0.0)
                    return "ERROR: An invalid scale factor was specified. The scale factor must be the decimal representation of a single-precision (32-bit) floating point number greater than zero.";
                iArgIndex++;
                break;

            case 'u':
            case 'U':
                iArgIndex++;
                *pUniformLines = 1;
                break;

            default:
                return "ERROR: Unrecognized command line argument.\n";
        }
    }

    // Return successfully.

    return NULL;
}

int main(int argc, char *argv[])
{
    printf("vectorshp version 1.1.0 (3-Feb-2012)\n");
    printf("\n");

    // Parse the command line arguments.

    if (argc == 1 || argc == 2 && argv[1][0] == '-' && argv[1][1] == '?' && argv[1][2] == 0)
    {
        PrintUsage();
        return 1;
    }

    char *pszXBinary = NULL;
    char *pszYBinary = NULL;
    char *pszOutFile = NULL;
    DataType iDataType = Int8_DataType;
    int iCols = 0;
    int iRows = 0;
    double dCellSize = 0.0;
    double dXLowerLeftCorner = 0.0;
    double dYLowerLeftCorner = 0.0;
    int bUseNoData = 0;
    int iNoDataValue = 0;
    double dNoDataValue = 0;
    double dScaleFactor = 1.0;
    int iUniformLines = 0;

    char *pszErrorMessage = ParseCommandLineArguments(argc, argv, &pszXBinary, &pszYBinary, &pszOutFile, &iDataType, &iCols, &iRows, &dCellSize, &dXLowerLeftCorner, &dYLowerLeftCorner, &bUseNoData, &iNoDataValue, &dNoDataValue, &dScaleFactor, &iUniformLines);
    if (pszErrorMessage != NULL)
    {
        printf(pszErrorMessage);
        printf("\n");
        printf("Type \"vectorshp -?\" for help.\n");
        return 1;
    }

    // Make sure the output shapefile does not exist already.

    size_t iLen = strlen(pszOutFile);

    if (iLen >= 4 &&
        pszOutFile[iLen - 4] == '.' && 
        (pszOutFile[iLen - 3] == 's' || pszOutFile[iLen - 3] == 'S') &&
        (pszOutFile[iLen - 2] == 'h' || pszOutFile[iLen - 2] == 'H') &&
        (pszOutFile[iLen - 1] == 'p' || pszOutFile[iLen - 1] == 'P'))
    {
        char *pszNewOutFile = (char *) calloc(iLen - 3, 1);
        if (pszNewOutFile == NULL)
        {
            printf("ERROR: Out of memory.\n");
            return 1;
        }
        pszOutFile = strncpy(pszNewOutFile, pszOutFile, iLen - 4);
    }

    iLen = strlen(pszOutFile);
    char *pszTemp = (char *) calloc(iLen + 5, 1);
    if (pszTemp == NULL)
    {
        printf("ERROR: Out of memory.\n");
        return 1;
    }
    strcpy(pszTemp, pszOutFile);
    pszTemp[iLen] = '.';

    pszTemp[iLen + 1] = 's';
    pszTemp[iLen + 2] = 'h';
    pszTemp[iLen + 3] = 'p';
    FILE *pTempFile = fopen(pszTemp, "r");
    if (pTempFile != NULL)
    {
        printf("ERROR: The file \"%s\" already exists. Please delete it or specify a different output file name.\n", pszTemp);
        return 1;
    }

    pszTemp[iLen + 1] = 's';
    pszTemp[iLen + 2] = 'h';
    pszTemp[iLen + 3] = 'x';
    pTempFile = fopen(pszTemp, "r");
    if (pTempFile != NULL)
    {
        printf("ERROR: The file \"%s\" already exists. Please delete it or specify a different output file name.\n", pszTemp);
        return 1;
    }

    pszTemp[iLen + 1] = 'd';
    pszTemp[iLen + 2] = 'b';
    pszTemp[iLen + 3] = 'f';
    pTempFile = fopen(pszTemp, "r");
    if (pTempFile != NULL)
    {
        printf("ERROR: The file \"%s\" already exists. Please delete it or specify a different output file name.\n", pszTemp);
        return 1;
    }

    // Open the binary files file for reading.

    HANDLE hXFile = CreateFile(pszXBinary, GENERIC_READ, FILE_SHARE_READ, 0, OPEN_EXISTING, FILE_ATTRIBUTE_NORMAL + FILE_FLAG_SEQUENTIAL_SCAN, 0);

    if (hXFile == INVALID_HANDLE_VALUE)
    {
        DWORD dwError = GetLastError();
        LPSTR lpMsgBuf = NULL;
        FormatMessage(FORMAT_MESSAGE_ALLOCATE_BUFFER | FORMAT_MESSAGE_FROM_SYSTEM | FORMAT_MESSAGE_IGNORE_INSERTS, NULL, dwError, MAKELANGID(LANG_NEUTRAL, SUBLANG_DEFAULT), (LPSTR) &lpMsgBuf, 0, NULL);
        if (lpMsgBuf != NULL)
            printf("ERROR: Failed to open file \"%s\" for reading. The operating system returned error %d: %s\n", pszXBinary, dwError, lpMsgBuf);
        else
            printf("ERROR: Failed to open file \"%s\" for reading.\n", pszXBinary);
        return 1;
    }

    printf("Opened \"%s\" for reading.\n", pszXBinary);

    HANDLE hYFile = CreateFile(pszYBinary, GENERIC_READ, FILE_SHARE_READ, 0, OPEN_EXISTING, FILE_ATTRIBUTE_NORMAL + FILE_FLAG_SEQUENTIAL_SCAN, 0);

    if (hYFile == INVALID_HANDLE_VALUE)
    {
        DWORD dwError = GetLastError();
        LPSTR lpMsgBuf = NULL;
        FormatMessage(FORMAT_MESSAGE_ALLOCATE_BUFFER | FORMAT_MESSAGE_FROM_SYSTEM | FORMAT_MESSAGE_IGNORE_INSERTS, NULL, dwError, MAKELANGID(LANG_NEUTRAL, SUBLANG_DEFAULT), (LPSTR) &lpMsgBuf, 0, NULL);
        if (lpMsgBuf != NULL)
            printf("ERROR: Failed to open file \"%s\" for reading. The operating system returned error %d: %s\n", pszYBinary, dwError, lpMsgBuf);
        else
            printf("ERROR: Failed to open file \"%s\" for reading.\n", pszYBinary);
        return 1;
    }

    printf("Opened \"%s\" for reading.\n", pszYBinary);

    // Verify that the binary files are the right sizes.

    LONGLONG liExpectedSize = (LONGLONG) iCols * (LONGLONG) iRows * (LONGLONG) DataTypeSizes[iDataType];
    LARGE_INTEGER liFileSize;

    if (!GetFileSizeEx(hXFile, &liFileSize))
    {
        DWORD dwError = GetLastError();
        LPSTR lpMsgBuf = NULL;
        FormatMessage(FORMAT_MESSAGE_ALLOCATE_BUFFER | FORMAT_MESSAGE_FROM_SYSTEM | FORMAT_MESSAGE_IGNORE_INSERTS, NULL, dwError, MAKELANGID(LANG_NEUTRAL, SUBLANG_DEFAULT), (LPSTR) &lpMsgBuf, 0, NULL);
        if (lpMsgBuf != NULL)
            printf("ERROR: Failed to retrieve the size of the file \"%s\". The operating system returned error %d: %s\n", pszXBinary, dwError, lpMsgBuf);
        else
            printf("ERROR: Failed to retrieve the size of the file \"%s\".\n", pszXBinary);
        return 1;
    }

    if (liFileSize.QuadPart != liExpectedSize)
    {
        printf("ERROR: The size of file \"%s\" is %I64i bytes, but it was expected to be %I64i bytes, given that it has %i columns, %i rows, and its data type requires %i bytes. Please check your parameters.\n", pszXBinary, liFileSize.QuadPart, liExpectedSize, iCols, iRows, DataTypeSizes[iDataType]);
        return 1;
    }

    if (!GetFileSizeEx(hYFile, &liFileSize))
    {
        DWORD dwError = GetLastError();
        LPSTR lpMsgBuf = NULL;
        FormatMessage(FORMAT_MESSAGE_ALLOCATE_BUFFER | FORMAT_MESSAGE_FROM_SYSTEM | FORMAT_MESSAGE_IGNORE_INSERTS, NULL, dwError, MAKELANGID(LANG_NEUTRAL, SUBLANG_DEFAULT), (LPSTR) &lpMsgBuf, 0, NULL);
        if (lpMsgBuf != NULL)
            printf("ERROR: Failed to retrieve the size of the file \"%s\". The operating system returned error %d: %s\n", pszYBinary, dwError, lpMsgBuf);
        else
            printf("ERROR: Failed to retrieve the size of the file \"%s\".\n", pszYBinary);
        return 1;
    }

    if (liFileSize.QuadPart != liExpectedSize)
    {
        printf("ERROR: The size of file \"%s\" is %I64i bytes, but it was expected to be %I64i bytes, given that it has %i columns, %i rows, and its data type requires %i bytes. Please check your parameters.\n", pszYBinary, liFileSize.QuadPart, liExpectedSize, iCols, iRows, DataTypeSizes[iDataType]);
        return 1;
    }

    // Create the output shapefile.

    printf("Writing lines to output shapefile \"%s\"...\n", pszOutFile);

    errno = 0;
    SHPHandle hOutShapeFile = NULL;

    hOutShapeFile = SHPCreate(pszOutFile, SHPT_ARC);

    if (hOutShapeFile == NULL)
    {
        if (errno != 0)
            printf("ERROR: Failed to create the output shapefile \"%s\". The shapelib.dll function SHPCreate() failed. On exit of that function, the global error number was set to %d: %s\n", pszOutFile, errno, strerror(errno));
        else
            printf("ERROR: Failed to create the output shapefile \"%s\". The shapelib.dll function SHPCreate() failed. Unfortunately that function does not return any detailed error information, and the global error number was set to zero (no error).\n", pszOutFile);
        return 1;
    }

    // Create the output DBF file.

    DBFHandle hOutDBFFile = DBFCreate(pszOutFile);
    if (hOutDBFFile == NULL)
    {
        if (errno != 0)
            printf("ERROR: Failed to create the output DBF file \"%s.dbf\". The shapelib.dll function DBFCreate() failed. On exit of that function, the global error number was set to %d: %s\n", pszOutFile, errno, strerror(errno));
        else
            printf("ERROR: Failed to create the output DBF file \"%s.dbf\". The shapelib.dll function DBFCreate() failed. Unfortunately that function does not return any detailed error information, and the global error number was set to zero (no error).\n", pszOutFile);
        SHPClose(hOutShapeFile);
        return 1;
    }

    DBFAddField(hOutDBFFile, "Magnitude", FTDouble, 19, 11);
    DBFAddField(hOutDBFFile, "Direction", FTDouble, 19, 11);

    // Read cells from the binary files and write lines to the shapefile.
    
    BYTE pXReadBuf[BUF_SIZE];
    BYTE pYReadBuf[BUF_SIZE];
    DWORD dwBytesToRead = sizeof(pXReadBuf);
    DWORD dwBytesRead = 0;
    DWORD dwReadIndex = 0;

    for (INT32 row = 0; row < iRows; row++)
    {
        for (INT32 col = 0; col < iCols; col++)
        {
            // If we exhausted the read buffers, refill them from the binary
            // files.

            if (dwReadIndex >= dwBytesRead)
            {
                if (!ReadFile(hXFile, pXReadBuf, dwBytesToRead, &dwBytesRead, NULL))
                {
                    DWORD dwError = GetLastError();
                    LPSTR lpMsgBuf = NULL;
                    FormatMessage(FORMAT_MESSAGE_ALLOCATE_BUFFER | FORMAT_MESSAGE_FROM_SYSTEM | FORMAT_MESSAGE_IGNORE_INSERTS, NULL, dwError, MAKELANGID(LANG_NEUTRAL, SUBLANG_DEFAULT), (LPSTR) &lpMsgBuf, 0, NULL);
                    if (lpMsgBuf != NULL)
                        printf("ERROR: Failed to read from the file \"%s\". The operating system returned error %d: %s\n", pszXBinary, dwError, lpMsgBuf);
                    else
                        printf("ERROR: Failed to read from the file \"%s\".\n", pszXBinary);
                    SHPClose(hOutShapeFile);
                    DBFClose(hOutDBFFile);
                    return 1;
                }
                if (dwBytesRead == 0)
                {
                    printf("ERROR: Reached the end of file \"%s\" prematurely. It is too small for the number of columns and rows and the data type specified.\n", pszXBinary);
                    SHPClose(hOutShapeFile);
                    DBFClose(hOutDBFFile);
                    return 1;
                }

                DWORD dwBytesRead2 = 0;
                if (!ReadFile(hYFile, pYReadBuf, dwBytesRead, &dwBytesRead2, NULL))
                {
                    DWORD dwError = GetLastError();
                    LPSTR lpMsgBuf = NULL;
                    FormatMessage(FORMAT_MESSAGE_ALLOCATE_BUFFER | FORMAT_MESSAGE_FROM_SYSTEM | FORMAT_MESSAGE_IGNORE_INSERTS, NULL, dwError, MAKELANGID(LANG_NEUTRAL, SUBLANG_DEFAULT), (LPSTR) &lpMsgBuf, 0, NULL);
                    if (lpMsgBuf != NULL)
                        printf("ERROR: Failed to read from the file \"%s\". The operating system returned error %d: %s\n", pszYBinary, dwError, lpMsgBuf);
                    else
                        printf("ERROR: Failed to read from the file \"%s\".\n", pszYBinary);
                    SHPClose(hOutShapeFile);
                    DBFClose(hOutDBFFile);
                    return 1;
                }
                if (dwBytesRead2 == 0)
                {
                    printf("ERROR: Reached the end of file \"%s\" prematurely. It is too small for the number of columns and rows and the data type specified.\n", pszYBinary);
                    SHPClose(hOutShapeFile);
                    DBFClose(hOutDBFFile);
                    return 1;
                }
                if (dwBytesRead2 != dwBytesRead)
                {
                    printf("ERROR: Reached the end of file \"%s\" prematurely. It is not the same size as file \"%s\".\n", pszYBinary, pszXBinary);
                    SHPClose(hOutShapeFile);
                    DBFClose(hOutDBFFile);
                    return 1;
                }

                dwReadIndex = 0;
            }

            // Only create a line if neither cell is NODATA.

            bool bNoData = false;

            if (bUseNoData)
                switch (iDataType)
                {
                    case Int8_DataType:
                        bNoData = ((INT32) *((INT8 *) (pXReadBuf + dwReadIndex))) == iNoDataValue || ((INT32) *((INT8 *) (pYReadBuf + dwReadIndex))) == iNoDataValue;
                        break;

                    case Uint8_DataType:
                        bNoData = ((INT32) *((UINT8 *) (pXReadBuf + dwReadIndex))) == iNoDataValue || ((INT32) *((UINT8 *) (pYReadBuf + dwReadIndex))) == iNoDataValue;
                        break;

                    case Int16_DataType:
                        bNoData = ((INT32) *((INT16 *) (pXReadBuf + dwReadIndex))) == iNoDataValue || ((INT32) *((INT16 *) (pYReadBuf + dwReadIndex))) == iNoDataValue;
                        break;

                    case Uint16_DataType:
                        bNoData = ((INT32) *((UINT16 *) (pXReadBuf + dwReadIndex))) == iNoDataValue || ((INT32) *((UINT16 *) (pYReadBuf + dwReadIndex))) == iNoDataValue;
                        break;

                    case Int32_DataType:
                        bNoData = *((INT32 *) (pXReadBuf + dwReadIndex)) == iNoDataValue || *((INT32 *) (pYReadBuf + dwReadIndex)) == iNoDataValue;
                        break;

                    case Float_DataType:
                        bNoData = ((DOUBLE) *((FLOAT *) (pXReadBuf + dwReadIndex))) == dNoDataValue || ((DOUBLE) *((FLOAT *) (pYReadBuf + dwReadIndex))) == dNoDataValue;
                        break;

                    case Double_DataType:
                        bNoData = *((DOUBLE *) (pXReadBuf + dwReadIndex)) == dNoDataValue || *((DOUBLE *) (pYReadBuf + dwReadIndex)) == dNoDataValue;
                        break;

                    default:
                        printf("ERROR: Programming error in this tool. Invalid DataType. Please contact the author of this tool for assistance.\n");
                        SHPClose(hOutShapeFile);
                        DBFClose(hOutDBFFile);
                        return 1;
                }

            if (!bNoData)
            {
                // Calculate the starting coordinates of the vector.

                double dXStart = dXLowerLeftCorner + dCellSize * (double) col + dCellSize / 2;
                double dYStart = dYLowerLeftCorner + dCellSize * (double) (iRows - row - 1) + dCellSize / 2;

                // Calculate the ending coordinates of the vector.

                double dXEnd = dXStart;
                double dYEnd = dYStart;

                switch (iDataType)
                {
                    case Int8_DataType:
                        dXEnd += ((DOUBLE) *((INT8 *) (pXReadBuf + dwReadIndex))) * dScaleFactor;
                        dYEnd += ((DOUBLE) *((INT8 *) (pYReadBuf + dwReadIndex))) * dScaleFactor;
                        break;

                    case Uint8_DataType:
                        dXEnd += ((DOUBLE) *((UINT8 *) (pXReadBuf + dwReadIndex))) * dScaleFactor;
                        dYEnd += ((DOUBLE) *((UINT8 *) (pYReadBuf + dwReadIndex))) * dScaleFactor;
                        break;

                    case Int16_DataType:
                        dXEnd += ((DOUBLE) *((INT16 *) (pXReadBuf + dwReadIndex))) * dScaleFactor;
                        dYEnd += ((DOUBLE) *((INT16 *) (pYReadBuf + dwReadIndex))) * dScaleFactor;
                        break;

                    case Uint16_DataType:
                        dXEnd += ((DOUBLE) *((UINT16 *) (pXReadBuf + dwReadIndex))) * dScaleFactor;
                        dYEnd += ((DOUBLE) *((UINT16 *) (pYReadBuf + dwReadIndex))) * dScaleFactor;
                        break;

                    case Int32_DataType:
                        dXEnd += ((DOUBLE) *((INT32 *) (pXReadBuf + dwReadIndex))) * dScaleFactor;
                        dYEnd += ((DOUBLE) *((INT32 *) (pYReadBuf + dwReadIndex))) * dScaleFactor;
                        break;

                    case Float_DataType:
                        dXEnd += ((DOUBLE) *((FLOAT *) (pXReadBuf + dwReadIndex))) * dScaleFactor;
                        dYEnd += ((DOUBLE) *((FLOAT *) (pYReadBuf + dwReadIndex))) * dScaleFactor;
                        break;

                    case Double_DataType:
                        dXEnd += *((DOUBLE *) (pXReadBuf + dwReadIndex)) * dScaleFactor;
                        dYEnd += *((DOUBLE *) (pYReadBuf + dwReadIndex)) * dScaleFactor;
                        break;

                    default:
                        printf("ERROR: Programming error in this tool. Invalid DataType. Please contact the author of this tool for assistance.\n");
                        SHPClose(hOutShapeFile);
                        DBFClose(hOutDBFFile);
                        return 1;
                }

                // Write the vector to the shapefile.

                double dMagnitude = sqrt(pow((dXEnd - dXStart) / dScaleFactor, 2) + pow((dYEnd - dYStart) / dScaleFactor, 2));
                double dDirection = atan2(dYEnd - dYStart, dXEnd - dXStart) / 3.1415926535897931 * 180.0;

                double fX[2] = {dXStart, dXEnd};
                double fY[2] = {dYStart, dYEnd};

                if (iUniformLines)
                {
                    double dLineLength = dCellSize * dScaleFactor;
                    fX[1] = dXStart + (dXEnd - dXStart) * (dLineLength / dMagnitude);
                    fY[1] = dYStart + (dYEnd - dYStart) * (dLineLength / dMagnitude);

                    //if (dXEnd >= dXStart)
                    //    fX[1] = dXStart + sqrt(pow(dMagnitude, 2) - pow(dYStart - dYEnd, 2)) * dCellSize * dScaleFactor;
                    //else
                    //    fX[1] = dXStart - sqrt(pow(dMagnitude, 2) - pow(dYStart - dYEnd, 2)) * dCellSize * dScaleFactor;

                    //if (dYEnd >= dYStart)
                    //    fY[1] = dYStart + sqrt(pow(dMagnitude, 2) - pow(dXStart - dXEnd, 2)) * dCellSize * dScaleFactor;
                    //else
                    //    fY[1] = dYStart - sqrt(pow(dMagnitude, 2) - pow(dXStart - dXEnd, 2)) * dCellSize * dScaleFactor;
                }

                SHPObject *shpObj = SHPCreateObject(SHPT_ARC, -1, 0, NULL, NULL, 2, &fX[0], &fY[0], NULL, NULL);

                int iShapeID = SHPWriteObject(hOutShapeFile, -1, shpObj);
                if (iShapeID == -1)
                {
                    if (errno != 0)
                        printf("ERROR: Failed to write line to output shapefile \"%s\". The shapelib.dll function SHPWriteObject() failed. On exit of that function, the global error number was set to %d: %s\n", pszOutFile, errno, strerror(errno));
                    else
                        printf("ERROR: Failed to write line to output shapefile \"%s\". The shapelib.dll function SHPWriteObject() failed. Unfortunately that function does not return any detailed error information, and the global error number was set to zero (no error).\n", pszOutFile);
                    SHPClose(hOutShapeFile);
                    DBFClose(hOutDBFFile);
                    return 1;
                }

                SHPDestroyObject(shpObj);

                // Write attributes to the DBF file.


                if (!DBFWriteDoubleAttribute(hOutDBFFile, iShapeID, 0, dMagnitude) ||
                    !DBFWriteDoubleAttribute(hOutDBFFile, iShapeID, 1, dDirection))
                {
                    if (errno != 0)
                        printf("ERROR: Failed to write an attribute to the output DBF file \"%s.dbf\". The shapelib.dll function DBFWriteDoubleAttribute() failed. On exit of that function, the global error number was set to %d: %s\n", pszOutFile, errno, strerror(errno));
                    else
                        printf("ERROR: Failed to write an attribute to the output DBF file \"%s.dbf\". The shapelib.dll function DBFWriteDoubleAttribute() failed. Unfortunately that function does not return any detailed error information, and the global error number was set to zero (no error).\n", pszOutFile);
                    SHPClose(hOutShapeFile);
                    DBFClose(hOutDBFFile);
                    return 1;
                }
            }

            // Go on to the next cell.

            switch (iDataType)
            {
                case Int8_DataType:
                case Uint8_DataType:
                    dwReadIndex++;
                    break;

                case Int16_DataType:
                case Uint16_DataType:
                    dwReadIndex += 2;
                    break;

                case Int32_DataType:
                case Float_DataType:
                    dwReadIndex += 4;
                    break;

                case Double_DataType:
                    dwReadIndex += 8;
                    break;

                default:
                    printf("ERROR: Programming error in this tool. Invalid DataType. Please contact the author of this tool for assistance.\n");
                    SHPClose(hOutShapeFile);
                    DBFClose(hOutDBFFile);
                    return 1;
            }
        }
    }

    // Close the files.

    printf("Done.\n");

    CloseHandle(hXFile);
    CloseHandle(hYFile);
    SHPClose(hOutShapeFile);
    DBFClose(hOutDBFFile);

    // Return successfully.

        return 0;
}