I created a plugin few years ago in Python. The plugin uses QGIS interpolation and works for QGIS 2.18. But because the QGIS interpolation creates raster (grid) and I need only Triangulation (TIN as lines) I checked the QGIS code and call function
int QgsTinInterpolator::interpolatePoint( double x, double y, double &result, QgsFeedback * ) via PyQGIS.
The Python code for 2.18 QGIS API:
ld_layer = QgsInterpolator.LayerData() ld_layer.vectorLayer = layer ld_layer.zCoordInterpolation = False ld_layer.interpolationAttribute = index ld_layer.mInputType = 1 ld2_layer = QgsInterpolator.LayerData() ld2_layer.vectorLayer = layer2 ld2_layer.zCoordInterpolation = False ld2_layer.interpolationAttribute = index ld2_layer.mInputType = 1 itp = QgsTINInterpolator(lds) itp.setExportTriangulationToFile(True) itp.setTriangulationFilePath(self.path) rect = lds.vectorLayer.extent() current_y_value = rect.yMaximum() - self.resolution / 2.0 current_x_value = rect.xMinimum() + self.resolution / 2.0 itp.interpolatePoint(current_x_value, current_y_value)
I created LayerData objects and then I called
itp.interpolatePoint(current_x_value, current_y_value) the Triangulation filed was created. It is because of the
void QgsTinInterpolator::initialize() function. This function is inside the qgstininterpolator.cpp and because it is private function I can't call that function via PyQGIS API, so I decided call
interpolatePoint() which calls
initialize() if initialization is not configured.
But now mainly because of in the plugin there is a performance need and because I am really interesting in exploring the world of C/C++ I decided few months ago learn C++ and rewrite the plugin. There are some changes in API:
- there is no
- there is no
- instead of there is
So I started hacking with QGIS and I wrote some code in C++:
QgsVectorLayer *exp5Layer = new QgsVectorLayer("/home/username/Downloads/data/testing_data/exp_q5_1.shp", "exp_q5", "ogr"); QgsVectorLayer *exp20Layer = new QgsVectorLayer("/home/username/Downloads/data/testing_data/exp_q20_1.shp", "exp_q20", "ogr"); const int index = 0; QgsInterpolator::LayerData ld5; QgsInterpolator::LayerData ld20; ld5.source = exp5Layer; ld5.interpolationAttribute = index; ld5.sourceType = QgsInterpolator::SourceType::SourceStructureLines; ld5.valueSource = QgsInterpolator::ValueSource::ValueZ; ld20.source = exp20Layer; ld20.interpolationAttribute = index; ld20.sourceType = QgsInterpolator::SourceType::SourceStructureLines; ld20.valueSource = QgsInterpolator::ValueSource::ValueZ; QList<QgsInterpolator::LayerData> lds; lds.append(ld5); lds.append(ld20); QgsTinInterpolator *itp = new QgsTinInterpolator(lds); QgsWkbTypes::Type wkbType = exp5Layer->wkbType(); QgsFields fields = exp5Layer->fields(); QgsVectorLayer *sink = new QgsVectorLayer("/home/username/Downloads/data/testing_data/triangulation.shp", "triangulation", "ogr"); itp->setTriangulationSink(sink); double current_y_value = rect.yMaximum() - 0.1 / 2.0; double current_x_value = rect.xMinimum() + 0.1 / 2.0; double interpolatedValue; QgsFeedback *feedback = new QgsFeedback(); itp->interpolatePoint(current_x_value, current_y_value, interpolatedValue, feedback); sink->~QgsVectorFileWriter();
The code is pretty similiar to Python code but of course there are changes in the API of QGIS 3.7. But the problem is the triangulation file (triangulation.shp) is created but there are not features.
It is because in
bool DualEdgeTriangulation::saveTriangulation( QgsFeatureSink *sink, QgsFeedback *feedback ) const function in DualEdgeInterpolation.cpp there is a cycle (at line 2997) and because of
mHalfEdge.size() == 0 addFeature function
sink->addFeature( edgeLineFeature, QgsFeatureSink::FastInsert ); (at line 3031) is not called).
There is no error and triangulation file is created but without features. If I run interpolation from the QGIS GUI with same parameters the triangulation file is created with features correctly.
Is there any problem with new API (like more changes that I didn't notice) or there is any other problem? I can't figure why
mHalfEdge is empty.
I continued with the investigating the problem and I found there was some big changes in
QgsTINInterpolator::initialize() function (commit). Maybe this can lead to solution or some clue.