QGIS API Documentation 3.30.0-'s-Hertogenbosch (f186b8efe0)
qgsmeshlayerprofilegenerator.cpp
Go to the documentation of this file.
1/***************************************************************************
2 qgsmeshlayerprofilegenerator.cpp
3 ---------------
4 begin : March 2022
5 copyright : (C) 2022 by Nyall Dawson
6 email : nyall dot dawson at gmail dot com
7 ***************************************************************************/
8
9/***************************************************************************
10 * *
11 * This program is free software; you can redistribute it and/or modify *
12 * it under the terms of the GNU General Public License as published by *
13 * the Free Software Foundation; either version 2 of the License, or *
14 * (at your option) any later version. *
15 * *
16 ***************************************************************************/
18#include "qgsprofilerequest.h"
19#include "qgscurve.h"
20#include "qgsmeshlayer.h"
22#include "qgsgeos.h"
23#include "qgsterrainprovider.h"
24#include "qgsmeshlayerutils.h"
25#include "qgslinesymbol.h"
26#include "qgsfillsymbol.h"
28#include "qgsprofilesnapping.h"
29#include "qgsprofilepoint.h"
30
31//
32// QgsMeshLayerProfileGenerator
33//
34
36{
37 return QStringLiteral( "mesh" );
38}
39
40QVector<QgsProfileIdentifyResults> QgsMeshLayerProfileResults::identify( const QgsProfilePoint &point, const QgsProfileIdentifyContext &context )
41{
42 const QVector<QgsProfileIdentifyResults> noLayerResults = QgsAbstractProfileSurfaceResults::identify( point, context );
43
44 // we have to make a new list, with the correct layer reference set
45 QVector<QgsProfileIdentifyResults> res;
46 res.reserve( noLayerResults.size() );
47 for ( const QgsProfileIdentifyResults &result : noLayerResults )
48 {
49 res.append( QgsProfileIdentifyResults( mLayer, result.results() ) );
50 }
51 return res;
52}
53
54//
55// QgsMeshLayerProfileGenerator
56//
57
59 : mId( layer->id() )
60 , mFeedback( std::make_unique< QgsFeedback >() )
61 , mProfileCurve( request.profileCurve() ? request.profileCurve()->clone() : nullptr )
62 , mSourceCrs( layer->crs() )
63 , mTargetCrs( request.crs() )
64 , mTransformContext( request.transformContext() )
65 , mOffset( layer->elevationProperties()->zOffset() )
66 , mScale( layer->elevationProperties()->zScale() )
67 , mLayer( layer )
68 , mStepDistance( request.stepDistance() )
69{
70 layer->updateTriangularMesh();
71 mTriangularMesh = *layer->triangularMesh();
72
73 mSymbology = qgis::down_cast< QgsMeshLayerElevationProperties * >( layer->elevationProperties() )->profileSymbology();
74 mLineSymbol.reset( qgis::down_cast< QgsMeshLayerElevationProperties * >( layer->elevationProperties() )->profileLineSymbol()->clone() );
75 mFillSymbol.reset( qgis::down_cast< QgsMeshLayerElevationProperties * >( layer->elevationProperties() )->profileFillSymbol()->clone() );
76}
77
79{
80 return mId;
81}
82
84
86{
87 if ( !mProfileCurve || mFeedback->isCanceled() )
88 return false;
89
90 // we need to transform the profile curve to the mesh's CRS
91 QgsGeometry transformedCurve( mProfileCurve->clone() );
92 mLayerToTargetTransform = QgsCoordinateTransform( mSourceCrs, mTargetCrs, mTransformContext );
93
94 try
95 {
96 transformedCurve.transform( mLayerToTargetTransform, Qgis::TransformDirection::Reverse );
97 }
98 catch ( QgsCsException & )
99 {
100 QgsDebugMsg( QStringLiteral( "Error transforming profile line to mesh CRS" ) );
101 return false;
102 }
103
104 if ( mFeedback->isCanceled() )
105 return false;
106
107 mResults = std::make_unique< QgsMeshLayerProfileResults >();
108 mResults->mLayer = mLayer;
109 mResults->copyPropertiesFromGenerator( this );
110
111 // we don't currently have any method to determine line->mesh intersection points, so for now we just sample at about 100(?) points over the line
112 const double curveLength = transformedCurve.length();
113
114 if ( !std::isnan( mStepDistance ) )
115 transformedCurve = transformedCurve.densifyByDistance( mStepDistance );
116 else
117 transformedCurve = transformedCurve.densifyByDistance( curveLength / 100 );
118
119 if ( mFeedback->isCanceled() )
120 return false;
121
122 for ( auto it = transformedCurve.vertices_begin(); it != transformedCurve.vertices_end(); ++it )
123 {
124 if ( mFeedback->isCanceled() )
125 return false;
126
127 QgsPoint point = ( *it );
128 const double height = heightAt( point.x(), point.y() );
129
130 try
131 {
132 point.transform( mLayerToTargetTransform );
133 }
134 catch ( QgsCsException & )
135 {
136 continue;
137 }
138 mResults->mRawPoints.append( QgsPoint( point.x(), point.y(), height ) );
139 }
140
141 if ( mFeedback->isCanceled() )
142 return false;
143
144 // convert x/y values back to distance/height values
145 QgsGeos originalCurveGeos( mProfileCurve.get() );
146 originalCurveGeos.prepareGeometry();
147 QString lastError;
148 for ( const QgsPoint &pixel : std::as_const( mResults->mRawPoints ) )
149 {
150 if ( mFeedback->isCanceled() )
151 return false;
152
153 const double distance = originalCurveGeos.lineLocatePoint( pixel, &lastError );
154
155 if ( !std::isnan( pixel.z() ) )
156 {
157 mResults->minZ = std::min( pixel.z(), mResults->minZ );
158 mResults->maxZ = std::max( pixel.z(), mResults->maxZ );
159 }
160 mResults->mDistanceToHeightMap.insert( distance, pixel.z() );
161 }
162
163 return true;
164}
165
167{
168 return mResults.release();
169}
170
172{
173 return mFeedback.get();
174}
175
176double QgsMeshLayerProfileGenerator::heightAt( double x, double y )
177{
178 return QgsMeshLayerUtils::interpolateZForPoint( mTriangularMesh, x, y ) * mScale + mOffset;
179}
180
Abstract base class for storage of elevation profiles.
QVector< QgsProfileIdentifyResults > identify(const QgsProfilePoint &point, const QgsProfileIdentifyContext &context) override
Identify results visible at the specified profile point.
Class for doing transforms between two map coordinate systems.
Custom exception class for Coordinate Reference System related exceptions.
Definition: qgsexception.h:66
Base class for feedback objects to be used for cancellation of something running in a worker thread.
Definition: qgsfeedback.h:45
A geometry is the spatial representation of a feature.
Definition: qgsgeometry.h:164
double length() const
Returns the planar, 2-dimensional length of geometry.
QgsGeometry densifyByDistance(double distance) const
Densifies the geometry by adding regularly placed extra nodes inside each segment so that the maximum...
Qgis::GeometryOperationResult transform(const QgsCoordinateTransform &ct, Qgis::TransformDirection direction=Qgis::TransformDirection::Forward, bool transformZ=false) SIP_THROW(QgsCsException)
Transforms this geometry as described by the coordinate transform ct.
QgsAbstractGeometry::vertex_iterator vertices_begin() const
Returns STL-style iterator pointing to the first vertex of the geometry.
QgsAbstractGeometry::vertex_iterator vertices_end() const
Returns STL-style iterator pointing to the imaginary vertex after the last vertex of the geometry.
Does vector analysis using the geos library and handles import, export, exception handling*.
Definition: qgsgeos.h:99
void prepareGeometry() override
Prepares the geometry, so that subsequent calls to spatial relation methods are much faster.
Definition: qgsgeos.cpp:252
double lineLocatePoint(const QgsPoint &point, QString *errorMsg=nullptr) const
Returns a distance representing the location along this linestring of the closest point on this lines...
Definition: qgsgeos.cpp:2814
~QgsMeshLayerProfileGenerator() override
QgsFeedback * feedback() const override
Access to feedback object of the generator (may be nullptr)
QgsMeshLayerProfileGenerator(QgsMeshLayer *layer, const QgsProfileRequest &request)
Constructor for QgsMeshLayerProfileGenerator.
QgsAbstractProfileResults * takeResults() override
Takes results from the generator.
bool generateProfile(const QgsProfileGenerationContext &context=QgsProfileGenerationContext()) override
Generate the profile (based on data stored in the class).
QString sourceId() const override
Returns a unique identifier representing the source of the profile.
QString type() const override
Returns the unique string identifier for the results type.
QVector< QgsProfileIdentifyResults > identify(const QgsProfilePoint &point, const QgsProfileIdentifyContext &context) override
Identify results visible at the specified profile point.
Represents a mesh layer supporting display of data on structured or unstructured meshes.
Definition: qgsmeshlayer.h:100
void updateTriangularMesh(const QgsCoordinateTransform &transform=QgsCoordinateTransform())
Gets native mesh and updates (creates if it doesn't exist) the base triangular mesh.
QgsMapLayerElevationProperties * elevationProperties() override
Returns the layer's elevation properties.
QgsTriangularMesh * triangularMesh(double minimumTriangleSize=0) const
Returns triangular mesh (nullptr before rendering or calling to updateMesh).
Point geometry type, with support for z-dimension and m-values.
Definition: qgspoint.h:49
void transform(const QgsCoordinateTransform &ct, Qgis::TransformDirection d=Qgis::TransformDirection::Forward, bool transformZ=false) override SIP_THROW(QgsCsException)
Transforms the geometry using a coordinate transform.
Definition: qgspoint.cpp:382
Q_GADGET double x
Definition: qgspoint.h:52
double y
Definition: qgspoint.h:53
Encapsulates the context in which an elevation profile is to be generated.
Encapsulates the context of identifying profile results.
Stores identify results generated by a QgsAbstractProfileResults object.
Encapsulates a point on a distance-elevation profile.
Encapsulates properties and constraints relating to fetching elevation profiles from different source...
#define QgsDebugMsg(str)
Definition: qgslogger.h:38
const QgsCoordinateReferenceSystem & crs