I'm using GDAL v2.1.3 with the SWIG bindings in my custom C# application.
I've created a method to do some calculation of my raster file:
public bool GdalCalculate(string input, string output, string formula,
double? minValue = null)
{
if (!File.Exists(input))
throw new FileNotFoundException("Can't find the input file", new
Exception("Working with " + input));
// 1. First copy the input file to create the output file
try
{
/*
-------------------------------------------------------------------- */
/* Get driver
*/
/*
-------------------------------------------------------------------- */
using (var drv = Gdal.GetDriverByName("GTiff"))
{
if (drv == null)
{
throw new Exception("Can't get GTiff driver");
}
/*
-------------------------------------------------------------------- */
/* Open dataset.
*/
/*
-------------------------------------------------------------------- */
using (var ds = Gdal.Open(input, Access.GA_ReadOnly))
{
if (ds == null)
{
throw new Exception("Can't open GDAL dataset: " +
input);
}
var options = new[] { "" };
using (var newDataset = drv.CreateCopy(output, ds, 0,
options, null, "Sample Data"))
{
if (newDataset == null)
{
throw new Exception("Can't create destination
dataset: " + output);
}
// Close input dataset:
ds.Dispose();
// 2. Loop through all pixels and perform formula on it:
using (var band = newDataset.GetRasterBand(1))
{
double noData;
int hasValue;
band.GetNoDataValue(out noData, out hasValue);
var sizeX = band.XSize;
var numLines = band.YSize;
var func = new Function("f(A) = " + formula);
// Loop through each line in turn.
for (var line = 0; line < numLines; line++)
{
var scanline = new float[sizeX];
// Read in data for the current line
var cplReturn = band.ReadRaster(0, line, sizeX,
1, scanline, sizeX, 1, 0, 0);
if (cplReturn != CPLErr.CE_None)
{
throw new Exception("band.ReadRaster
failed: " + Gdal.GetLastErrorMsg());
}
var outputLine = new List<float>();
foreach (var f in scanline)
{
var pixelValue = f;
if ((float)f != (float)noData)
{
// Calculate
pixelValue = (float)func.calculate(f);
if (minValue.HasValue)
pixelValue =
(float)Math.Max(pixelValue, minValue.GetValueOrDefault());
}
outputLine.Add(pixelValue);
}
// Rewrite line:
cplReturn = band.WriteRaster(0, line, sizeX, 1,
outputLine.ToArray(), sizeX, 1, 0, 0);
if (cplReturn != CPLErr.CE_None)
{
throw new Exception("band.WriteRaster
failed: " + Gdal.GetLastErrorMsg());
}
}
// 3. Save changes to file:
band.FlushCache();
newDataset.FlushCache();
}
}
}
}
}
catch (Exception e)
{
throw new Exception("Can't open GDAL dataset: " + input, e);
}
return true;
}
This method is working fine, but is very slow.
I'm using a 3621x4466 tiff file, which has a size of 67MB.
It is taking around 10 minutes to process.
(float)func.calculate(f); is using http://mathparser.org/
An example of a function I'm using: 9.7E-05-e^(3.1 + 0.9 * A) where A is
the pixel value.
How can I optimize my function?
Thanks,
Paul
_______________________________________________
gdal-dev mailing list
gdal-dev@lists.osgeo.org
https://lists.osgeo.org/mailman/listinfo/gdal-dev
