GDAL
GNM Tutorial

This document is intended to describe using the GNM C++ classes to work with networks. It is advised to read the GNM architecture before to understand the purpose and structure of GNM classes.

Managing networks

In the first example we will create a small water network on the base of the set of spatial data (two shapefiles: pipes and wells which are situated at the GDAL source tree: autotest\gnm\data). The use of the common network format - GNMGdalNetwork class - will allow us to select one of the GDAL-supported vector formats for our network - ESRI Shapefile. After the creation we will build a topology and add some additional data: pumps layer, in order to manually edit network topology.

Initially we register GDAL drivers and create some options (string pairs), which will be passed as parameters during network creation. Here we create a network's name.

#include "gnm.h"
#include <vector>
int main ()
{
char **papszDSCO = NULL;
papszDSCO = CSLAddNameValue(papszDSCO, GNM_MD_NAME, "my_pipes_network");
papszDSCO = CSLAddNameValue(papszDSCO, GNM_MD_SRS, "EPSG:4326");
papszDSCO = CSLAddNameValue(papszDSCO, GNM_MD_DESCR, "My pipes network");
papszDSCO = CSLAddNameValue(papszDSCO, GNM_MD_FORMAT, "ESRI Shapefile");

Some options are obligatory. The following parameters must be specified during the network creation: the path/name; format of network storage; spatial reference system (EPSG, WKT, etc.). The according dataset with the "network part" will be created and the resulting network will be returned.

GNMGenericNetwork* poDS = (GNMGenericNetwork*) poDriver->Create( "..\\network_data", 0, 0, 0, GDT_Unknown,
papszDSCO );
CSLDestroy(papszDSCO);

For now we have a void network consisted of only "system layers". We need to populate it with "class layers" full of features, so we open a certain foreign dataset and copy layers from it to our network. Note, that we use GDALDataset:: methods for working with "class layers", because GNMNetwork inherited from GDALDataset.

GDALDataset *poSrcDS = (GDALDataset*) GDALOpenEx("..\\in_data",
GDAL_OF_VECTOR | GDAL_OF_READONLY, NULL, NULL, NULL );
OGRLayer *poSrcLayer1 = poSrcDS->GetLayerByName("pipes");
OGRLayer *poSrcLayer2 = poSrcDS->GetLayerByName("wells");
poDS->CopyLayer(poSrcLayer1, "pipes");
poDS->CopyLayer(poSrcLayer2, "wells");
GDALClose(poSrcDS);

After the successful copying we have the network full of features, but with no topology. The features were added and registered in the network but they are still not connected with each other. Now it is time to build the network topology. There are two ways of doing this in GNM: manually or automatically. In the most cases automatic building is more convenient, while manual is useful for small editings. Automatic building requires some parameters: we must specify which "class layers" will participate in topology building (we select our two layers), a snap tolerance, direct and inverse cost, direction, which is equal 0.00005 in our case. If the building will be successful the network's graph will be filled with the according connections.

printf("\nBuilding network topology ...\n");
char **papszLayers = NULL;
for(int i = 0; i < poDS->GetLayerCount(); ++i)
{
OGRLayer* poLayer = poDS->GetLayer(i);
papszLayers = CSLAddString(papszLayers, poLayer->GetName() );
}
if(poGenericNetwork->ConnectPointsByLines(papszLayers, dfTolerance,
dfDirCost, dfInvCost, eDir) != CE_None )
{
printf("Building topology failed\n");
}
else
{
printf("Topology has been built successfully\n");
}

At this point we have a ready network with topological and spatial data, which can be used now for different purposes (analysis, converting into different formats, etc). But sometimes it is necessary to modify some network's data. For example we need to add additional features and attach them to our built topology (modify topology). We create a new "class layer" in the network and add one feature to it.

OGRLayer *poNewLayer = poDS->CreateLayer("pumps", , NULL, wkbPoint, NULL );
if( poNewLayer == NULL )
{
printf( "Layer creation failed.\n" );
exit( 1 );
}
OGRFieldDefn fieldDefn ("pressure",OFTReal);
if( poNewLayer->CreateField( &fieldDefn ) != OGRERR_NONE )
{
printf( "Creating Name field failed.\n" );
exit( 1 );
}
pt.setX(37.291466);
pt.setY(55.828351);
poFeature->SetGeometry(&pt);
if( poNewLayer->CreateFeature( poFeature ) != OGRERR_NONE )
{
printf( "Failed to create feature.\n" );
exit( 1 );
}
GNMGFID gfid = poFeature->GetFID();

After the successful creation the feature will be registered in the network and we can connect it with others. There can be two possible ways to do this. In the first case we need a real feature which will be an edge in the connection, while in the second case we do not need such feature, and passing -1 into the ConnectFeatures() method means that the special system edge will be created for this connection and added to the graph automatically. In our case we had added only one point feature and we have not got the line one to be an edge, so we will use the "virtual" connection. We pass the GFID of our point as the source, the GFID of one of the existed features as the target and -1 as the connector. Note that we also set the costs (direct and inverse) and the direction of our edge manually and these values will be written to the graph. When we used the automatic connection (which also uses ConnectFeatures() internally) such vales were set automatically according to the rule which we also set before.

if (poDS->ConnectFeatures(gfid ,63, -1, 5.0, 5.0, GNMDirection_SrcToTgt) != GNMError_None)
{
printf("Can not connect features\n");
}

After all we correctly close the network which frees the allocated resources.

GDALClose(poDS);

All in one block:

#include "gnm.h"
#include "gnm_priv.h"
int main ()
{
char **papszDSCO = NULL;
papszDSCO = CSLAddNameValue(papszDSCO, GNM_MD_NAME, "my_pipes_network");
papszDSCO = CSLAddNameValue(papszDSCO, GNM_MD_SRS, "EPSG:4326");
papszDSCO = CSLAddNameValue(papszDSCO, GNM_MD_DESCR, "My pipes network");
papszDSCO = CSLAddNameValue(papszDSCO, GNM_MD_FORMAT, "ESRI Shapefile");
GDALDriver *poDriver = GetGDALDriverManager()->GetDriverByName("GNMFile");
GNMGenericNetwork* poDS = (GNMGenericNetwork*) poDriver->Create( "..\\network_data", 0, 0, 0, GDT_Unknown,
papszDSCO );
CSLDestroy(papszDSCO);
if (poDS == NULL)
{
printf("Failed to create network\n");
exit(1);
}
GDALDataset *poSrcDS = (GDALDataset*) GDALOpenEx("..\\in_data",GDAL_OF_VECTOR | GDAL_OF_READONLY, NULL, NULL, NULL );
if(poSrcDS == NULL)
{
printf("Can not open source dataset at\n");
exit(1);
}
OGRLayer *poSrcLayer1 = poSrcDS->GetLayerByName("pipes");
OGRLayer *poSrcLayer2 = poSrcDS->GetLayerByName("wells");
if (poSrcLayer1 == NULL || poSrcLayer2 == NULL)
{
printf("Can not process layers of source dataset\n");
exit(1);
}
poDS->CopyLayer(poSrcLayer1, "pipes");
poDS->CopyLayer(poSrcLayer2, "wells");
GDALClose(poSrcDS);
printf("\nBuilding network topology ...\n");
char **papszLayers = NULL;
for(int i = 0; i < poDS->GetLayerCount(); ++i)
{
OGRLayer* poLayer = poDS->GetLayer(i);
papszLayers = CSLAddString(papszLayers, poLayer->GetName() );
}
if(poGenericNetwork->ConnectPointsByLines(papszLayers, dfTolerance,
dfDirCost, dfInvCost, eDir) != CE_None )
{
printf("Building topology failed\n");
exit(1);
}
else
{
printf("Topology has been built successfully\n");
}
OGRLayer *poNewLayer = poDS->CreateLayer("pumps", , NULL, wkbPoint, NULL );
if( poNewLayer == NULL )
{
printf( "Layer creation failed.\n" );
exit( 1 );
}
OGRFieldDefn fieldDefn ("pressure",OFTReal);
if( poNewLayer->CreateField( &fieldDefn ) != OGRERR_NONE )
{
printf( "Creating Name field failed.\n" );
exit( 1 );
}
OGRFeature *poFeature = OGRFeature::CreateFeature(poNewLayer->GetLayerDefn());
pt.setX(37.291466);
pt.setY(55.828351);
poFeature->SetGeometry(&pt);
if( poNewLayer->CreateFeature( poFeature ) != OGRERR_NONE )
{
printf( "Failed to create feature.\n" );
exit( 1 );
}
GNMGFID gfid = poFeature->GetFID();
if (poDS->ConnectFeatures(gfid ,63, -1, 5.0, 5.0, GNMDirection_SrcToTgt) != GNMError_None)
{
printf("Can not connect features\n");
}
GDALClose(poDS);
}

Analysing networks

In the second example we will analyse the network which we have built in the first example. We will calculate the shortest path between two points via Dijkstra algorithm performing the feature blockings and saving the resulting path into the file.

Initially we open our network, passing the path to its Shapefile dataset.

#include "gnm.h"
#include "gnm_priv.h"
int main ()
{
GNMGenericNetwork *poNet = (GNMGenericNetwork*) GDALOpenEx("..\\network_data",GDAL_OF_GNM | GDAL_OF_UPDATE, NULL, NULL, NULL );
if(poSrcDS == NULL)
{
printf("Can not open source dataset at\n");
exit(1);
}

Before any calculations we open the dataset which will hold the layer with the resulting path.

GDALDataset *poResDS;
poResDS = (GDALDataset*) GDALOpenEx("..\\out_data",
NULL, NULL, NULL);
if (poResDS == NULL)
{
printf("Failed to open resulting dataset\n");
exit(1);
}

Finally we use the Dijkstra shortest path method to calculations. This path will be found passing over the blocked feature and saved into internal memory OGRLayer, which we copy to the real dataset. Now it can be visualized by GIS.

OGRLayer *poResLayer = poNet->GetPath(64, 41, GATDijkstraShortestPath, NULL);
if (poResLayer == NULL)
{
printf("Failed to save or calculate path\n");
}
else if (poResDS->CopyLayer(poResLayer, "shp_tutorial.shp") == NULL)
{
printf("Failed to save path to the layer\n");
}
else
{
printf("Path saved successfully\n");
}
GDALClose(poResDS);
poNet->ReleaseResultSet(poRout);
GDALClose(poNet);
}

All in one block:

#include "gnm.h"
#include "gnmstdanalysis.h"
int main ()
{
GNMGenericNetwork *poNet = (GNMGenericNetwork*) GDALOpenEx("..\\network_data",
NULL, NULL, NULL );
if(poSrcDS == NULL)
{
printf("Can not open source dataset at\n");
exit(1);
}
GDALDataset *poResDS;
poResDS = (GDALDataset*) GDALOpenEx("..\\out_data",
NULL, NULL, NULL);
if (poResDS == NULL)
{
printf("Failed to open resulting dataset\n");
exit(1);
}
poNet->ChangeBlockState(36, true);
OGRLayer *poResLayer = poNet->GetPath(64, 41, GATDijkstraShortestPath, NULL);
if (poResLayer == NULL)
{
printf("Failed to save or calculate path\n");
}
else if (poResDS->CopyLayer(poResLayer, "shp_tutorial.shp") == NULL)
{
printf("Failed to save path to the layer\n");
}
else
{
printf("Path saved successfully\n");
}
GDALClose(poResDS);
poNet->ReleaseResultSet(poRout);
GDALClose(poNet);
}
OGRLayer::GetName
virtual const char * GetName()
Return the layer name.
Definition: ogrlayer.cpp:1727
OGRLayer::CreateField
virtual OGRErr CreateField(OGRFieldDefn *poField, int bApproxOK=TRUE)
Create a new field on a layer.
Definition: ogrlayer.cpp:664
GNMGenericNetwork::GetLayerCount
virtual int GetLayerCount() override
Get the number of layers in this dataset.
GDALDataset::CreateLayer
virtual OGRLayer * CreateLayer(const char *pszName, OGRSpatialReference *poSpatialRef=nullptr, OGRwkbGeometryType eGType=wkbUnknown, char **papszOptions=nullptr)
This method attempts to create a new layer on the dataset with the indicated name,...
Definition: gdaldataset.cpp:4341
CSLAddString
char ** CSLAddString(char **papszStrList, const char *pszNewString) CPL_WARN_UNUSED_RESULT
Append a string to a StringList and return a pointer to the modified StringList.
Definition: cpl_string.cpp:81
OGRFeature::SetGeometry
OGRErr SetGeometry(const OGRGeometry *)
Set feature geometry.
Definition: ogrfeature.cpp:436
wkbPoint
0-dimensional geometric object, standard WKB
Definition: ogr_core.h:320
GDT_Unknown
Definition: gdal.h:60
GDALDriver
Format specific driver.
Definition: gdal_priv.h:1422
GDALClose
void GDALClose(GDALDatasetH)
Close GDAL dataset.
Definition: gdaldataset.cpp:3588
OGRLayer
This class represents a layer of simple features, with access methods.
Definition: ogrsf_frmts.h:69
OGRPoint
Point class.
Definition: ogr_geometry.h:809
GDALDataset::CopyLayer
virtual OGRLayer * CopyLayer(OGRLayer *poSrcLayer, const char *pszNewName, char **papszOptions=nullptr)
Duplicate an existing layer.
Definition: gdaldataset.cpp:4786
OGRLayer::GetLayerDefn
virtual OGRFeatureDefn * GetLayerDefn()=0
Fetch the schema information for this layer.
GDALDataset
A set of associated raster bands, usually from one file.
Definition: gdal_priv.h:334
OGRPoint::setX
void setX(double xIn)
Set x.
Definition: ogr_geometry.h:865
GNMGenericNetwork::CopyLayer
virtual OGRLayer * CopyLayer(OGRLayer *poSrcLayer, const char *pszNewName, char **papszOptions=nullptr) override
Duplicate an existing layer.
GDALDriverManager::GetDriverByName
GDALDriver * GetDriverByName(const char *)
Fetch a driver based on the short name.
Definition: gdaldrivermanager.cpp:589
GDALDataset::ReleaseResultSet
virtual void ReleaseResultSet(OGRLayer *poResultsSet)
Release results of ExecuteSQL().
Definition: gdaldataset.cpp:6470
GNMGenericNetwork::ChangeBlockState
virtual CPLErr ChangeBlockState(GNMGFID nFID, bool bIsBlock)
Change the block state of edge or vertex.
OGRLayer::CreateFeature
OGRErr CreateFeature(OGRFeature *poFeature) CPL_WARN_UNUSED_RESULT
Create and write a new feature within a layer.
Definition: ogrlayer.cpp:625
GDAL_OF_READONLY
#define GDAL_OF_READONLY
Open in read-only mode.
Definition: gdal.h:423
GDALOpenEx
GDALDatasetH GDALOpenEx(const char *pszFilename, unsigned int nOpenFlags, const char *const *papszAllowedDrivers, const char *const *papszOpenOptions, const char *const *papszSiblingFiles) CPL_WARN_UNUSED_RESULT
Open a raster or vector file as a GDALDataset.
Definition: gdaldataset.cpp:3206
OGRFieldDefn
Definition of an attribute of an OGRFeatureDefn.
Definition: ogr_feature.h:91
OGRFeature::CreateFeature
static OGRFeature * CreateFeature(OGRFeatureDefn *)
Feature factory.
Definition: ogrfeature.cpp:245
GDALDataset::GetLayerByName
virtual OGRLayer * GetLayerByName(const char *)
Fetch a layer by name.
Definition: gdaldataset.cpp:5179
GNMGenericNetwork::GetPath
virtual OGRLayer * GetPath(GNMGFID nStartFID, GNMGFID nEndFID, GNMGraphAlgorithmType eAlgorithm, char **papszOptions) override
Create path between start and end GFIDs.
GNMGenericNetwork::ConnectFeatures
virtual CPLErr ConnectFeatures(GNMGFID nSrcFID, GNMGFID nTgtFID, GNMGFID nConFID=-1, double dfCost=1, double dfInvCost=1, GNMDirection eDir=GNM_EDGE_DIR_BOTH)
Connects two features via third feature (may be virtual, so the identificator should be -1).
GDALDriver::Create
GDALDataset * Create(const char *pszName, int nXSize, int nYSize, int nBands, GDALDataType eType, char **papszOptions) CPL_WARN_UNUSED_RESULT
Create a new dataset with this driver.
Definition: gdaldriver.cpp:161
GDAL_OF_UPDATE
#define GDAL_OF_UPDATE
Open in update mode.
Definition: gdal.h:429
GDALAllRegister
void GDALAllRegister(void)
Register all known configured GDAL drivers.
Definition: gdalallregister.cpp:61
CSLDestroy
void CSLDestroy(char **papszStrList)
Free string list.
Definition: cpl_string.cpp:198
OGRFeature
A simple feature, including geometry and attributes.
Definition: ogr_feature.h:353
CSLAddNameValue
char ** CSLAddNameValue(char **papszStrList, const char *pszName, const char *pszValue) CPL_WARN_UNUSED_RESULT
Add a new entry to a StringList of "Name=Value" pairs, ("Name:Value" pairs are also supported for bac...
Definition: cpl_string.cpp:1822
OGRERR_NONE
#define OGRERR_NONE
Success.
Definition: ogr_core.h:291
OGRFeature::DestroyFeature
static void DestroyFeature(OGRFeature *)
Destroy feature.
Definition: ogrfeature.cpp:281
GetGDALDriverManager
GDALDriverManager * GetGDALDriverManager(void)
Fetch the global GDAL driver manager.
Definition: gdaldrivermanager.cpp:96
OFTReal
Double Precision floating point.
Definition: ogr_core.h:597
GDAL_OF_VECTOR
#define GDAL_OF_VECTOR
Allow vector drivers to be used.
Definition: gdal.h:447
GNMGenericNetwork::GetLayer
virtual OGRLayer * GetLayer(int) override
Fetch a layer by index.
GNMGenericNetwork
GNM class which represents a geography network of generic format.
Definition: gnm.h:194
OGRFeature::GetFID
GIntBig GetFID() const
Get feature identifier.
Definition: ogr_feature.h:711
GDAL_OF_GNM
#define GDAL_OF_GNM
Allow gnm drivers to be used.
Definition: gdal.h:453

Generated for GDAL by doxygen 1.8.16.