How to translate Tiff.ReadEncodedTile to elevation

Question

I'm new with working with reading tiff images and I'm trying to get the elevation terrain values from a tiff map by using LibTiff. The maps I need to decode are tile organized. Below the fragment of the code I'm using currently to get these values, based on the library documentation and research on the web:

    private void getBytes()
    {
        int numBytes = bitsPerSample / 8;           //Number of bytes depending the tiff map
        int stride = numBytes * height;
        byte[] bufferTiff = new byte[stride * height];  // this is the buffer with the tiles data

        int offset = 0;

        for (int i = 0; i < tif.NumberOfTiles() - 1; i++)
        {
            int rawTileSize = (int)tif.RawTileSize(i);
            offset += tif.ReadEncodedTile(i, bufferTiff, offset, rawTileSize); 

        }

        values = new double[height, width];         // this is the matrix to save the heigth values in meters

        int ptr = 0;                    // pointer for saving the each data bytes
        int m = 0;
        int n = 0;

        byte[] byteValues = new byte[numBytes];     // bytes of each height data

        for (int i = 0; i < bufferTiff.Length; i++)
        {
            byteValues[ptr] = bufferTiff[i];

            ptr++;
            if (ptr % numBytes == 0)
            {
                ptr = 0;

                    if (n == height) // tiff map Y pixels
                    {
                        n = 0;
                        m++;
                        if (m == width) // tiff map X pixels
                        {
                            m = 0;
                        }
                    }

                    values[m, n] = BitConverter.ToDouble(byteValues, 0);    // Converts each byte data to the height value in meters. If the map is 32 bps the method I use is BitConverter.ToFloat

                    if (n == height - 1 && m == width - 1)
                        break;
                    n++;

            }
        }
        SaveArrayAsCSV(values, "values.txt");               
    }

    //Only to show results in a cvs file:
    public void SaveArrayAsCSV(double[,] arrayToSave, string fileName)  // source: http://stackoverflow.com/questions/8666518/how-can-i-write-a-general-array-to-csv-file
    {
        using (StreamWriter file = new StreamWriter(fileName))
        {
            WriteItemsToFile(arrayToSave, file);
        }
    }

   //Only to show results in a cvs file:
    private void WriteItemsToFile(Array items, TextWriter file)     // source: http://stackoverflow.com/questions/8666518/how-can-i-write-a-general-array-to-csv-file
    {
        int cont = 0;
        foreach (object item in items)
        {
            if (item is Array)
            {
                WriteItemsToFile(item as Array, file);
                file.Write(Environment.NewLine);
            }
            else {
                file.Write(item + " | ");
                cont++;
                if(cont == width)                       
                {
                    file.Write("\n");
                    cont = 0;
                }
            }
        }
    }

I've been testing two different maps (32 and 64 bits per sample) and the results are similar: At the begining, the data seems to be consistent, but there is a point in which all the other values are corrupted (even zero at the end of data results). I deduce there are some bytes that need to be ignored, but I don't know how identify them to depurate my code. The method Tiff.ReadScanline does not work for me, because the maps I need to decode are tiles organized, and this method is not for working with these kind of images (according to BitMiracle.LibTiff documentation). The method Tiff.ReadRGBATile is not valid neither, because the tiff images are not RGB. I can read these values with Matlab, but my project needs to be built in C#, so I can compare the expected results with mine. As reference (I think it could be helpful), these are some data extracted from one of the tiff files with LibTiff tag reading methods:

• ImageWidth: 2001
• ImageLength: 2001
• BitsPerSample: 32
• Compression: PackBits (aka Macintosh RLE)
• Photometric: MinIsBlack
• SamplesPerPixel: 1
• PlanarConfig: Contig
• TileWidth: 208
• TileLength: 208
• SampleFormat: 3

Thanks in advance by your help guys!

Answer 1

Here is an improved code, works with non-square tiles:

int imageWidth = tiff.GetField(TiffTag.IMAGEWIDTH)[0].ToInt();
                int imageHeight = tiff.GetField(TiffTag.IMAGELENGTH)[0].ToInt();
                int bytesPerSample = (int)tiff.GetField(TiffTag.BITSPERSAMPLE)[0].ToInt() / 8;
                SampleFormat format = (SampleFormat)tiff.GetField(TiffTag.SAMPLEFORMAT)[0].ToInt();

                //Array to return
                float[,] decoded = new float[imageHeight, imageWidth];

                //Get decode function (I only want a float array)
                Func<byte[], int, float> decode = GetConversionFunction(format, bytesPerSample);
                if (decode == null)
                {
                    throw new ArgumentException("Unsupported TIFF format:"+format);
                }

                if(tiff.IsTiled())
                {
                    //tile dimensions in pixels - the image dimensions MAY NOT be a multiple of these dimensions
                    int tileWidth = tiff.GetField(TiffTag.TILEWIDTH)[0].ToInt();
                    int tileHeight = tiff.GetField(TiffTag.TILELENGTH)[0].ToInt();

                    //tile matrix size
                    int numTiles = tiff.NumberOfTiles();
                    int tileMatrixWidth = (int)Math.Ceiling(imageWidth / (float)tileWidth);
                    int tileMatrixHeight = (int)Math.Ceiling(imageHeight / (float)tileHeight);

                    //tile dimensions in bytes
                    int tileBytesWidth = tileWidth * bytesPerSample;
                    int tileBytesHeight = tileHeight * bytesPerSample;

                    //tile buffer
                    int tileBufferSize = tiff.TileSize();
                    byte[] tileBuffer = new byte[tileBufferSize];

                    int imageHeightMinus1 = imageHeight - 1;

                    for (int tileIndex = 0 ; tileIndex < numTiles; tileIndex++)
                    {
                        int tileX = tileIndex / tileMatrixWidth;
                        int tileY = tileIndex % tileMatrixHeight;

                        tiff.ReadTile(tileBuffer, 0, tileX*tileWidth, tileY*tileHeight, 0, 0);

                        int xImageOffset = tileX * tileWidth;
                        int yImageOffset = tileY * tileHeight;

                        for (int col = 0; col < tileWidth && xImageOffset+col < imageWidth; col++ )
                        {
                            for(int row = 0; row < tileHeight && yImageOffset+row < imageHeight; row++)
                            {
                                decoded[imageHeightMinus1-(yImageOffset+row), xImageOffset+col] = decode(tileBuffer, row * tileBytesWidth + col * bytesPerSample);

                            }
                        }
                    }
                }

Answer 2

Ok, Finally I found the solution: My mistake was the parameter "count" in the function Tiff.ReadEncodedTile(tile, buffer, offset, count). The Tiff.RawTileSize(int) function, returns the compressed bytes size of the tile (different for each tile, depending of the compression algorithm), but Tiff.ReadEncodedTile returns the decompressed bytes (bigger and constant for all tiles). That's why not all the information was been saved properly, but just a part of data. Below the correct code with the terrain elevation matrix (need optimization but it works, I think it could be helpful)

 private void getBytes()
    {
        int numBytes = bitsPerSample / 8;
        int numTiles = tif.NumberOfTiles();
        int stride = numBytes * height;
        int bufferSize = tileWidth * tileHeight * numBytes * numTiles;
        int bytesSavedPerTile = tileWidth * tileHeight * numBytes; //this is the real size of the decompressed bytes
        byte[] bufferTiff = new byte[bufferSize];

        FieldValue[] value = tif.GetField(TiffTag.TILEWIDTH);
        int tilewidth = value[0].ToInt();

        value = tif.GetField(TiffTag.TILELENGTH);
        int tileHeigth = value[0].ToInt();

        int matrixSide = (int)Math.Sqrt(numTiles); // this works for a square image (for example a tiles organized tiff image)
        int bytesWidth = matrixSide * tilewidth;
        int bytesHeigth = matrixSide * tileHeigth;

        int offset = 0;

        for (int j = 0; j < numTiles; j++)
        {
            offset += tif.ReadEncodedTile(j, bufferTiff, offset, bytesSavedPerTile); //Here was the mistake. Now it works!
        }

        double[,] aux = new double[bytesHeigth, bytesWidth]; //Double for a 64 bps tiff image. This matrix will save the alldata, including the transparency (the "blank zone" I was talking before)

        terrainElevation = new double[height, width]; // Double for a 64 bps tiff image. This matrix will save only the elevation values, without transparency

        int ptr = 0;
        int m = 0;
        int n = -1;
        int contNumTile = 1;
        int contBytesPerTile = 0;
        int i = 0;
        int tileHeigthReference = tileHeigth;
        int tileWidthReference = tileWidth;
        int row = 1;
        int col = 1;

        byte[] bytesHeigthMeters = new byte[numBytes]; // Buffer to save each one elevation value to parse

        while (i < bufferTiff.Length && contNumTile < numTiles + 1)
        {
            for (contBytesPerTile = 0; contBytesPerTile < bytesSavedPerTile; contBytesPerTile++)
            {
                bytesHeigthMeters[ptr] = bufferTiff[i];
                ptr++;
                if (ptr % numBytes == 0 && ptr != 0)
                {
                    ptr = 0;
                    n++;

                    if (n == tileHeigthReference)
                    {
                        n = tileHeigthReference - tileHeigth;
                        m++;
                        if (m == tileWidthReference)
                        {
                            m = tileWidthReference - tileWidth;
                        }
                    }
                    double heigthMeters = BitConverter.ToDouble(bytesHeigthMeters, 0);

                    if (n < bytesWidth)
                    {
                        aux[m, n] = heigthMeters;
                    }
                    else
                    {
                        n = -1;
                    }

                }
                i++;
            }

            if (i % tilewidth == 0)
            {
                col++;
                if (col == matrixSide + 1)
                {
                    col = 1;
                }
            }

            if (contNumTile % matrixSide == 0)
            {
                row++;
                n = -1;
                if (row == matrixSide + 1)
                {
                    row = 1;

                }
            }

            contNumTile++;
            tileHeigthReference = tileHeight * (col);
            tileWidthReference = tileWidth * (row);

            m = tileWidth * (row - 1);

        }

        for (int x = 0; x < height; x++)
        {
            for (int y = 0; y < width; y++)
            {
                terrainElevation[x, y] = aux[x, y]; // Final result. Each position of matrix has saved each pixel terrain elevation of the map
            }
        }

    }

Regards!