I am trying to use libtiff.net to read elevation data from a GeoTIFF file. So far I have mostly just been able to read metadata from the file using the example at libtiff.net's webpage.
But howto read elevation data I do not understand... I tried first reading with Tiff.ReadScanline() as described here but the file I have seems to be stored differently (probably in tiles if I understand it correctly)
Here is the metadata (as far as I have been able to read) (the tiff file is from the danish terrain elevation data set):
Tiff c:\Users***\DTM_1km_6170_500.tif, page 0 has following tags set:
IMAGEWIDTH System.Int32 : 2500
IMAGELENGTH System.Int32 : 2500
BITSPERSAMPLE System.Int16 : 32
COMPRESSION BitMiracle.LibTiff.Classic.Compression : ADOBE_DEFLATE
PHOTOMETRIC BitMiracle.LibTiff.Classic.Photometric : MINISBLACK
STRIPOFFSETS System.UInt64[] : System.UInt64[]
SAMPLESPERPIXEL System.Int16 : 1
STRIPBYTECOUNTS System.UInt64[] : System.UInt64[]
PLANARCONFIG BitMiracle.LibTiff.Classic.PlanarConfig : CONTIG
PREDICTOR BitMiracle.LibTiff.Classic.Predictor : FLOATINGPOINT
TILEWIDTH System.Int32 : 256
TILELENGTH System.Int32 : 256
TILEOFFSETS System.UInt64[] : System.UInt64[]
TILEBYTECOUNTS System.UInt64[] : System.UInt64[]
SAMPLEFORMAT BitMiracle.LibTiff.Classic.SampleFormat : IEEEFP
DATATYPE System.Int16 : 3
GEOTIFF_MODELPIXELSCALETAG System.Int32 : 3 GEOTIFF_MODELPIXELSCALETAG System.Byte[] : Ù?Ù?
GEOTIFF_MODELTIEPOINTTAG System.Int32 : 6 GEOTIFF_MODELTIEPOINTTAG System.Byte[] : A ^WA
34735 System.Int32 : 36 34735 System.Byte[] : ± ± #° èd )#
34736 System.Int32 : 3 34736 System.Byte[] :
34737 System.Int32 : 30 34737 System.Byte[] : ETRS89 / UTM zone 32N|ETRS89|
42113 System.Int32 : 6 42113 System.Byte[] : -9999
The code I have written so far is as follows:
namespace GeoTIFFReader
{
public class GeoTIFF
{
private double[,] heightmap;
private double dx;
private double dy;
private double startx;
private double starty;
public GeoTIFF(string fn)
{
using (Tiff tiff = Tiff.Open(fn, "r"))
{
if (tiff == null)
{
// Error - could not open
return;
}
int width = tiff.GetField(TiffTag.IMAGEWIDTH)[0].ToInt();
int height = tiff.GetField(TiffTag.IMAGELENGTH)[0].ToInt();
heightmap = new double[width, height];
FieldValue[] modelPixelScaleTag = tiff.GetField(TiffTag.GEOTIFF_MODELPIXELSCALETAG);
FieldValue[] modelTiePointTag = tiff.GetField(TiffTag.GEOTIFF_MODELTIEPOINTTAG);
byte[] modelPixelScale = modelPixelScaleTag[1].GetBytes();
dx = BitConverter.ToDouble(modelPixelScale, 0);
dy = BitConverter.ToDouble(modelPixelScale, 8) * -1;
byte[] modelTransformation = modelTiePointTag[1].GetBytes();
double originLon = BitConverter.ToDouble(modelTransformation, 24);
double originLat = BitConverter.ToDouble(modelTransformation, 32);
startx = originLon + dx / 2.0;
starty = originLat + dy / 2.0;
double curx = startx;
double cury = starty;
FieldValue[] bitsPerSampleTag = tiff.GetField(TiffTag.BITSPERSAMPLE);
FieldValue[] tilewtag = tiff.GetField(TiffTag.TILEWIDTH);
FieldValue[] tilehtag = tiff.GetField(TiffTag.TILELENGTH);
int tilew = tilewtag[0].ToInt();
int tileh = tilehtag[0].ToInt();
var tile = new byte[tilew*tileh];
//var scanline = new byte[tiff.ScanlineSize()]; Does not work... wrong format
for (int il = 0; il < height; il++)
{
//tiff.ReadScanline(scanline, il); // Load il'th line of data
for (int ir = 0; ir < width; ir++)
{
// Here I would like to read each pixel data that contains elevation in gray-scale in f32 as far I as I understand from metadata
//object value = scanline[ir];
//heightmap[ir, il] = double.Parse(value.ToString());
}
}
Console.WriteLine(heightmap.ToString());
}
}
}
}
So if anyone knows howto extract this data, that would be very much appreciated.
So I stumbled across some hinting that lead me to find an answer to the specific question..:
int tileSize = tiff.TileSize();
for (int iw = 0; iw < nWidth; iw += tilew)
{
for (int ih = 0; ih < nHeight; ih += tileh)
{
byte[] buffer = new byte[tileSize];
tiff.ReadTile(buffer, 0, iw, ih, 0, 0);
for (int itw = 0; itw < tilew; itw++)
{
int iwhm = ih + itw;
if (iwhm > nWidth - 1)
{
break;
}
for (int ith = 0; ith < tileh; ith++)
{
int iyhm = iw + ith;
if (iyhm > nHeight - 1)
{
break;
}
heightMap[iwhm, iyhm] =
BitConverter.ToSingle(buffer, (itw * tileh + ith) * 4);
}
}
}
}
EDIT 2018-09-20:
@Graviton - Sorry for the long response time.. but here is the complete class I have used with a constructor that takes the filename as input and populates the heightMap (but @Nazonokaizijin looks a bit nicer and slimmer):
using System;
using BitMiracle.LibTiff.Classic;
using System.IO;
namespace GeoTIFFReader
{
public class GeoTIFF
{
private float[,] heightMap;
public float[,] HeightMap
{
get { return heightMap; }
private set { heightMap = value; }
}
private int nWidth;
public int NWidth
{
get { return nWidth; }
private set { nWidth = value; }
}
private int nHeight;
public int NHeight
{
get { return nHeight; }
private set { nHeight = value; }
}
private double dW;
public double DW
{
get { return dW; }
private set { dW = value; }
}
private double dH;
public double DH
{
get { return dH; }
private set { dH = value; }
}
private double startW;
public double StartW
{
get { return startW; }
private set { startW = value; }
}
private double startH;
public double StartH
{
get { return startH; }
private set { startH = value; }
}
public GeoTIFF(string fn)
{
using (Tiff tiff = Tiff.Open(fn, "r"))
{
if (tiff == null)
{
// Error - could not open
return;
}
nWidth = tiff.GetField(TiffTag.IMAGEWIDTH)[0].ToInt();
nHeight = tiff.GetField(TiffTag.IMAGELENGTH)[0].ToInt();
heightMap = new float[nWidth, nHeight];
FieldValue[] modelPixelScaleTag = tiff.GetField(TiffTag.GEOTIFF_MODELPIXELSCALETAG);
FieldValue[] modelTiePointTag = tiff.GetField(TiffTag.GEOTIFF_MODELTIEPOINTTAG);
byte[] modelPixelScale = modelPixelScaleTag[1].GetBytes();
dW = BitConverter.ToDouble(modelPixelScale, 0);
dH = BitConverter.ToDouble(modelPixelScale, 8) * -1;
byte[] modelTransformation = modelTiePointTag[1].GetBytes();
double originLon = BitConverter.ToDouble(modelTransformation, 24);
double originLat = BitConverter.ToDouble(modelTransformation, 32);
startW = originLon + dW / 2.0;
startH = originLat + dH / 2.0;
FieldValue[] tileByteCountsTag = tiff.GetField(TiffTag.TILEBYTECOUNTS);
long[] tileByteCounts = tileByteCountsTag[0].TolongArray();
FieldValue[] bitsPerSampleTag = tiff.GetField(TiffTag.BITSPERSAMPLE);
int bytesPerSample = bitsPerSampleTag[0].ToInt() / 8;
FieldValue[] tilewtag = tiff.GetField(TiffTag.TILEWIDTH);
FieldValue[] tilehtag = tiff.GetField(TiffTag.TILELENGTH);
int tilew = tilewtag[0].ToInt();
int tileh = tilehtag[0].ToInt();
int tileWidthCount = nWidth / tilew;
int remainingWidth = nWidth - tileWidthCount * tilew;
if (remainingWidth > 0)
{
tileWidthCount++;
}
int tileHeightCount = nHeight / tileh;
int remainingHeight = nHeight - tileHeightCount * tileh;
if (remainingHeight > 0)
{
tileHeightCount++;
}
int tileSize = tiff.TileSize();
for (int iw = 0; iw < nWidth; iw += tilew)
{
for (int ih = 0; ih < nHeight; ih += tileh)
{
byte[] buffer = new byte[tileSize];
tiff.ReadTile(buffer, 0, iw, ih, 0, 0);
for (int itw = 0; itw < tilew; itw++)
{
int iwhm = ih + itw;
if (iwhm > nWidth - 1)
{
break;
}
for (int ith = 0; ith < tileh; ith++)
{
int iyhm = iw + ith;
if (iyhm > nHeight - 1)
{
break;
}
heightMap[iwhm, iyhm] =
BitConverter.ToSingle(buffer, (itw * tileh + ith) * 4);
}
}
}
}
}
}
}
}