coordinatesdemelevationqgisqgis-2
I am working in QGIS 2.6. I have a file of several polygons and I want for each one, to find the lowest point (where is it located and its elevation).
I am using "Grid statistics for polygons" which gives me the elevation but not the coordinates of this lowest point.
How can I find these coordinates? Is there an existing tool in QGIS to perform this?
On the bottom of the application you can set extent or the center of the map.
It looks like this:
l
The first button switches between extent and map center.
Thanks to Vince's comments I was able to perform succesfully the whole process.
As I wrote, first I wanted to establish a number of AOIs on the raster. Nothing fancy here, I got raster properties by using this small function below:
public static IRasterProps GetRasterProperties(IRasterDataset rasterDataset, int rasterBandIndex)
{
IRasterBandCollection rasterBands = (IRasterBandCollection)rasterDataset;
var rasterBand = rasterBands.Item(rasterBandIndex);
return (IRasterProps)rasterBand;
}
Given the raster properties, I created a number of areas each defined by an IEnvelope. In order to split the input raster to AOIs, I used the ExtractByRectangle Geoprocessing tool.
ExtractByRectangle extract = new ExtractByRectangle(inputRaster, envelope, path);
IGeoProcessorResult2 result = gp.Execute(extract, null) as IGeoProcessorResult2;
parameters:
inputRaster - IRasterDataset (input raster)
envelope - IEnvelope (definition of our AOI)
path - file path to resulting raster
I won't describe all Geoprocessing shenanigans, but that is where ArcObjects get very cryptic in my opinion. A number of error codes that have no obvious explanation doesn't help too. I found great code by Kirk Kuykendall that helps a lot with debugging here: Avoiding fails from ArcObjects geoprocessing with .NET?
Now, we have our small raster AOIs. I need to get the point of maximal elevation for each of them. Hence, I start with computing statistics for every AOI using the function below:
public void ComputeRasterStatistics(IRasterDataset rasterDataset, int rasterBandIndex)
{
IRasterBandCollection rasterBands = (IRasterBandCollection)rasterDataset;
var rasterBand = rasterBands.Item(rasterBandIndex);
RasterStatistics = rasterBand.Statistics;
}
Statistics deliver information concerning values of raster extremes, inclusing the maximum point value of a raster. Next I convert the raster to points using the function below.
public static IFeatureClass RasterToPoints(IRasterDataset raster)
{
IConversionOp conversionOp = new RasterConversionOpClass();
IWorkspace shapeWS = FeatureWorkspaceHelper.CreateInMemoryWorkspace();
var featClass = conversionOp.RasterDataToPointFeatureData((IGeoDataset)raster, shapeWS, Guid.NewGuid().ToString());
return (IFeatureClass)featClass;
}
FeatureWorkspaceHelper.CreateInMemoryWorkspace(); is a helper function of mine, that creates an empty in-memory workspace. Have in mind this can be obviously done with 'classic' file workspace too.
Now I just need to find a point that has the maximal elevation and return it. The function below (some hardcodes there!) does that.
private IPoint GetRasterMaxElevationPoint(IFeatureClass featureClass, double val, int elevationIndex)
{
IQueryFilter queryFilter=new QueryFilterClass();
queryFilter.WhereClause = "GRID_CODE >= " + (val - 0.01).ToString();
var cursor = featureClass.Search(queryFilter, true);
IFeature feature = null;
IGeometry shape = null;
double maxValue = double.MinValue;
while ((feature = cursor.NextFeature()) != null)
{
if ((double) feature.Value[elevationIndex] > maxValue)
{
shape = feature.Shape;
maxValue = (double) feature.Value[elevationIndex];
}
}
return new PointClass()
{
X = shape.Envelope.LowerLeft.X,
Y = shape.Envelope.LowerLeft.Y,
SpatialReference = shape.SpatialReference
};
}
Best Answer
I can propose a (lengthy) solution using existing tools in QGIS.
Assumption:
Steps:
[Warning: This output point file can become very large for large areas or high resolution elevation data; and will make your computer non-responsive. Please try on a smaller area / dataset first.]
The output from step 2 should have geometry columns (lat & long), ID and minimum elevation of the corresponding polygons. A simple attribute query or sorting based on elevation will help you select the point having minimum elevation in each polygon.