QGIS API Documentation 3.40.0-Bratislava (b56115d8743)
Loading...
Searching...
No Matches
qgspointcloudlayerchunkloader_p.cpp
Go to the documentation of this file.
1/***************************************************************************
2 qgspointcloudlayerchunkloader_p.cpp
3 --------------------------------------
4 Date : October 2020
5 Copyright : (C) 2020 by Peter Petrik
6 Email : zilolv dot sk at gmail dot com
7 ***************************************************************************
8 * *
9 * This program is free software; you can redistribute it and/or modify *
10 * it under the terms of the GNU General Public License as published by *
11 * the Free Software Foundation; either version 2 of the License, or *
12 * (at your option) any later version. *
13 * *
14 ***************************************************************************/
15
17
18#include "qgs3dutils.h"
20#include "qgschunknode.h"
21#include "qgslogger.h"
22#include "qgspointcloudindex.h"
24#include "qgseventtracing.h"
25
26
31
32#include <QtConcurrent>
33#if QT_VERSION < QT_VERSION_CHECK(6, 0, 0)
34#include <Qt3DRender/QAttribute>
35#else
36#include <Qt3DCore/QAttribute>
37#endif
38#include <Qt3DRender/QTechnique>
39#include <Qt3DRender/QShaderProgram>
40#include <Qt3DRender/QGraphicsApiFilter>
41#include <QPointSize>
42
44
45
47
48QgsPointCloudLayerChunkLoader::QgsPointCloudLayerChunkLoader( const QgsPointCloudLayerChunkLoaderFactory *factory, QgsChunkNode *node, std::unique_ptr< QgsPointCloud3DSymbol > symbol,
49 const QgsCoordinateTransform &coordinateTransform, double zValueScale, double zValueOffset )
50 : QgsChunkLoader( node )
51 , mFactory( factory )
52 , mContext( factory->mRenderContext, coordinateTransform, std::move( symbol ), zValueScale, zValueOffset )
53{
54
55 QgsPointCloudIndex *pc = mFactory->mPointCloudIndex;
56 mContext.setAttributes( pc->attributes() );
57
58 const QgsChunkNodeId nodeId = node->tileId();
59 const IndexedPointCloudNode pcNode( nodeId.d, nodeId.x, nodeId.y, nodeId.z );
60
61 Q_ASSERT( pc->hasNode( pcNode ) );
62
63 QgsDebugMsgLevel( QStringLiteral( "loading entity %1" ).arg( node->tileId().text() ), 2 );
64
65 if ( mContext.symbol()->symbolType() == QLatin1String( "single-color" ) )
66 mHandler.reset( new QgsSingleColorPointCloud3DSymbolHandler() );
67 else if ( mContext.symbol()->symbolType() == QLatin1String( "color-ramp" ) )
68 mHandler.reset( new QgsColorRampPointCloud3DSymbolHandler() );
69 else if ( mContext.symbol()->symbolType() == QLatin1String( "rgb" ) )
70 mHandler.reset( new QgsRGBPointCloud3DSymbolHandler() );
71 else if ( mContext.symbol()->symbolType() == QLatin1String( "classification" ) )
72 {
73 mHandler.reset( new QgsClassificationPointCloud3DSymbolHandler() );
74 const QgsClassificationPointCloud3DSymbol *classificationSymbol = dynamic_cast<const QgsClassificationPointCloud3DSymbol *>( mContext.symbol() );
75 mContext.setFilteredOutCategories( classificationSymbol->getFilteredOutCategories() );
76 }
77
78 //
79 // this will be run in a background thread
80 //
81 mFutureWatcher = new QFutureWatcher<void>( this );
82 connect( mFutureWatcher, &QFutureWatcher<void>::finished, this, &QgsChunkQueueJob::finished );
83
84 const QgsAABB bbox = node->bbox();
85 const QFuture<void> future = QtConcurrent::run( [pc, pcNode, bbox, this]
86 {
87 const QgsEventTracing::ScopedEvent e( QStringLiteral( "3D" ), QStringLiteral( "PC chunk load" ) );
88
89 if ( mContext.isCanceled() )
90 {
91 QgsDebugMsgLevel( QStringLiteral( "canceled" ), 2 );
92 return;
93 }
94
95 mHandler->processNode( pc, pcNode, mContext );
96
97 if ( mContext.isCanceled() )
98 {
99 QgsDebugMsgLevel( QStringLiteral( "canceled" ), 2 );
100 return;
101 }
102
103 if ( mContext.symbol()->renderAsTriangles() )
104 mHandler->triangulate( pc, pcNode, mContext, bbox );
105 } );
106
107 // emit finished() as soon as the handler is populated with features
108 mFutureWatcher->setFuture( future );
109}
110
111QgsPointCloudLayerChunkLoader::~QgsPointCloudLayerChunkLoader()
112{
113 if ( mFutureWatcher && !mFutureWatcher->isFinished() )
114 {
115 disconnect( mFutureWatcher, &QFutureWatcher<void>::finished, this, &QgsChunkQueueJob::finished );
116 mContext.cancelRendering();
117 mFutureWatcher->waitForFinished();
118 }
119}
120
121void QgsPointCloudLayerChunkLoader::cancel()
122{
123 mContext.cancelRendering();
124}
125
126Qt3DCore::QEntity *QgsPointCloudLayerChunkLoader::createEntity( Qt3DCore::QEntity *parent )
127{
128 QgsPointCloudIndex *pc = mFactory->mPointCloudIndex;
129 const QgsChunkNodeId nodeId = mNode->tileId();
130 const IndexedPointCloudNode pcNode( nodeId.d, nodeId.x, nodeId.y, nodeId.z );
131 Q_ASSERT( pc->hasNode( pcNode ) );
132
133 Qt3DCore::QEntity *entity = new Qt3DCore::QEntity( parent );
134 mHandler->finalize( entity, mContext );
135 return entity;
136}
137
139
140
141QgsPointCloudLayerChunkLoaderFactory::QgsPointCloudLayerChunkLoaderFactory( const Qgs3DRenderContext &context, const QgsCoordinateTransform &coordinateTransform, QgsPointCloudIndex *pc, QgsPointCloud3DSymbol *symbol,
142 double zValueScale, double zValueOffset, int pointBudget )
143 : mRenderContext( context )
144 , mCoordinateTransform( coordinateTransform )
145 , mPointCloudIndex( pc )
146 , mZValueScale( zValueScale )
147 , mZValueOffset( zValueOffset )
148 , mPointBudget( pointBudget )
149{
150 mSymbol.reset( symbol );
151
152 try
153 {
154 mExtent = mCoordinateTransform.transformBoundingBox( mRenderContext.extent(), Qgis::TransformDirection::Reverse );
155 }
156 catch ( const QgsCsException & )
157 {
158 // bad luck, can't reproject for some reason
159 QgsDebugError( QStringLiteral( "Transformation of extent failed." ) );
160 }
161}
162
163QgsChunkLoader *QgsPointCloudLayerChunkLoaderFactory::createChunkLoader( QgsChunkNode *node ) const
164{
165 const QgsChunkNodeId id = node->tileId();
166
167 Q_ASSERT( mPointCloudIndex->hasNode( IndexedPointCloudNode( id.d, id.x, id.y, id.z ) ) );
168 QgsPointCloud3DSymbol *symbol = static_cast< QgsPointCloud3DSymbol * >( mSymbol->clone() );
169 return new QgsPointCloudLayerChunkLoader( this, node, std::unique_ptr< QgsPointCloud3DSymbol >( symbol ), mCoordinateTransform, mZValueScale, mZValueOffset );
170}
171
172int QgsPointCloudLayerChunkLoaderFactory::primitivesCount( QgsChunkNode *node ) const
173{
174 const QgsChunkNodeId id = node->tileId();
175 const IndexedPointCloudNode n( id.d, id.x, id.y, id.z );
176 Q_ASSERT( mPointCloudIndex->hasNode( n ) );
177 return mPointCloudIndex->nodePointCount( n );
178}
179
180QgsAABB nodeBoundsToAABB( QgsPointCloudDataBounds nodeBounds, QgsVector3D offset, QgsVector3D scale, const Qgs3DRenderContext &context, const QgsCoordinateTransform &coordinateTransform, double zValueOffset );
181
182QgsChunkNode *QgsPointCloudLayerChunkLoaderFactory::createRootNode() const
183{
184 const QgsAABB bbox = nodeBoundsToAABB( mPointCloudIndex->nodeBounds( IndexedPointCloudNode( 0, 0, 0, 0 ) ), mPointCloudIndex->offset(), mPointCloudIndex->scale(), mRenderContext, mCoordinateTransform, mZValueOffset );
185 const float error = mPointCloudIndex->nodeError( IndexedPointCloudNode( 0, 0, 0, 0 ) );
186 QgsChunkNode *node = new QgsChunkNode( QgsChunkNodeId( 0, 0, 0, 0 ), bbox, error );
187 node->setRefinementProcess( mSymbol->renderAsTriangles() ? Qgis::TileRefinementProcess::Replacement : Qgis::TileRefinementProcess::Additive );
188 return node;
189}
190
191QVector<QgsChunkNode *> QgsPointCloudLayerChunkLoaderFactory::createChildren( QgsChunkNode *node ) const
192{
193 QVector<QgsChunkNode *> children;
194 const QgsChunkNodeId nodeId = node->tileId();
195 const QgsAABB bbox = node->bbox();
196 const float childError = node->error() / 2;
197 float xc = bbox.xCenter(), yc = bbox.yCenter(), zc = bbox.zCenter();
198
199 for ( int i = 0; i < 8; ++i )
200 {
201 int dx = i & 1, dy = !!( i & 2 ), dz = !!( i & 4 );
202 const QgsChunkNodeId childId( nodeId.d + 1, nodeId.x * 2 + dx, nodeId.y * 2 + dy, nodeId.z * 2 + dz );
203
204 if ( !mPointCloudIndex->hasNode( IndexedPointCloudNode( childId.d, childId.x, childId.y, childId.z ) ) )
205 continue;
206 if ( !mExtent.isEmpty() &&
207 !mPointCloudIndex->nodeMapExtent( IndexedPointCloudNode( childId.d, childId.x, childId.y, childId.z ) ).intersects( mExtent ) )
208 continue;
209
210 // the Y and Z coordinates below are intentionally flipped, because
211 // in chunk node IDs the X,Y axes define horizontal plane,
212 // while in our 3D scene the X,Z axes define the horizontal plane
213 const float chXMin = dx ? xc : bbox.xMin;
214 const float chXMax = dx ? bbox.xMax : xc;
215 // Z axis: values are increasing to the south
216 const float chZMin = !dy ? zc : bbox.zMin;
217 const float chZMax = !dy ? bbox.zMax : zc;
218 const float chYMin = dz ? yc : bbox.yMin;
219 const float chYMax = dz ? bbox.yMax : yc;
220 QgsChunkNode *child = new QgsChunkNode( childId, QgsAABB( chXMin, chYMin, chZMin, chXMax, chYMax, chZMax ), childError, node );
221 child->setRefinementProcess( mSymbol->renderAsTriangles() ? Qgis::TileRefinementProcess::Replacement : Qgis::TileRefinementProcess::Additive );
222 children << child;
223 }
224 return children;
225}
226
228
229
230QgsAABB nodeBoundsToAABB( QgsPointCloudDataBounds nodeBounds, QgsVector3D offset, QgsVector3D scale, const Qgs3DRenderContext &context, const QgsCoordinateTransform &coordinateTransform, double zValueOffset )
231{
232 QgsVector3D extentMin3D( nodeBounds.xMin() * scale.x() + offset.x(), nodeBounds.yMin() * scale.y() + offset.y(), nodeBounds.zMin() * scale.z() + offset.z() + zValueOffset );
233 QgsVector3D extentMax3D( nodeBounds.xMax() * scale.x() + offset.x(), nodeBounds.yMax() * scale.y() + offset.y(), nodeBounds.zMax() * scale.z() + offset.z() + zValueOffset );
234 QgsCoordinateTransform extentTransform = coordinateTransform;
235 extentTransform.setBallparkTransformsAreAppropriate( true );
236 try
237 {
238 extentMin3D = extentTransform.transform( extentMin3D );
239 extentMax3D = extentTransform.transform( extentMax3D );
240 }
241 catch ( QgsCsException & )
242 {
243 QgsDebugError( QStringLiteral( "Error transforming node bounds coordinate" ) );
244 }
245 const QgsVector3D worldExtentMin3D = Qgs3DUtils::mapToWorldCoordinates( extentMin3D, context.origin() );
246 const QgsVector3D worldExtentMax3D = Qgs3DUtils::mapToWorldCoordinates( extentMax3D, context.origin() );
247 QgsAABB rootBbox( worldExtentMin3D.x(), worldExtentMin3D.y(), worldExtentMin3D.z(),
248 worldExtentMax3D.x(), worldExtentMax3D.y(), worldExtentMax3D.z() );
249 return rootBbox;
250}
251
252
253QgsPointCloudLayerChunkedEntity::QgsPointCloudLayerChunkedEntity( Qgs3DMapSettings *map, QgsPointCloudIndex *pc,
254 const QgsCoordinateTransform &coordinateTransform, QgsPointCloud3DSymbol *symbol,
255 float maximumScreenSpaceError, bool showBoundingBoxes,
256 double zValueScale, double zValueOffset,
257 int pointBudget )
258 : QgsChunkedEntity( map,
259 maximumScreenSpaceError,
260 new QgsPointCloudLayerChunkLoaderFactory( Qgs3DRenderContext::fromMapSettings( map ), coordinateTransform, pc, symbol, zValueScale, zValueOffset, pointBudget ), true, pointBudget )
261{
262 setShowBoundingBoxes( showBoundingBoxes );
263}
264
265QgsPointCloudLayerChunkedEntity::~QgsPointCloudLayerChunkedEntity()
266{
267 // cancel / wait for jobs
268 cancelActiveJobs();
269}
270
271QVector<QgsRayCastingUtils::RayHit> QgsPointCloudLayerChunkedEntity::rayIntersection( const QgsRayCastingUtils::Ray3D &ray, const QgsRayCastingUtils::RayCastContext &context ) const
272{
273 QVector<QgsRayCastingUtils::RayHit> result;
274 QgsPointCloudLayerChunkLoaderFactory *factory = static_cast<QgsPointCloudLayerChunkLoaderFactory *>( mChunkLoaderFactory );
275
276 // transform ray
277 const QgsVector3D rayOriginMapCoords = factory->mRenderContext.worldToMapCoordinates( ray.origin() );
278 const QgsVector3D pointMapCoords = factory->mRenderContext.worldToMapCoordinates( ray.origin() + ray.origin().length() * ray.direction().normalized() );
279 QgsVector3D rayDirectionMapCoords = pointMapCoords - rayOriginMapCoords;
280 rayDirectionMapCoords.normalize();
281
282 const int screenSizePx = std::max( context.screenSize.height(), context.screenSize.width() );
283
284 const QgsPointCloud3DSymbol *symbol = factory->mSymbol.get();
285 // Symbol can be null in case of no rendering enabled
286 if ( !symbol )
287 return result;
288 const double pointSize = symbol->pointSize();
289
290 // We're using the angle as a tolerance, effectively meaning we're fetching points intersecting a cone.
291 // This may be revisited to use a cylinder instead, if the balance between near/far points does not scale
292 // well with different point sizes, screen sizes and fov values.
293 const double limitAngle = 2. * pointSize / screenSizePx * factory->mRenderContext.fieldOfView();
294
295 // adjust ray to elevation properties
296 const QgsVector3D adjustedRayOrigin = QgsVector3D( rayOriginMapCoords.x(), rayOriginMapCoords.y(), ( rayOriginMapCoords.z() - factory->mZValueOffset ) / factory->mZValueScale );
297 QgsVector3D adjustedRayDirection = QgsVector3D( rayDirectionMapCoords.x(), rayDirectionMapCoords.y(), rayDirectionMapCoords.z() / factory->mZValueScale );
298 adjustedRayDirection.normalize();
299
300 QgsPointCloudIndex *index = factory->mPointCloudIndex;
301
302 const QgsPointCloudAttributeCollection attributeCollection = index->attributes();
303 QgsPointCloudRequest request;
304 request.setAttributes( attributeCollection );
305
306 double minDist = -1.;
307 const QList<QgsChunkNode *> activeNodes = this->activeNodes();
308 for ( QgsChunkNode *node : activeNodes )
309 {
310 const QgsChunkNodeId id = node->tileId();
311 const IndexedPointCloudNode n( id.d, id.x, id.y, id.z );
312
313 if ( !index->hasNode( n ) )
314 continue;
315
316 if ( !QgsRayCastingUtils::rayBoxIntersection( ray, node->bbox() ) )
317 continue;
318
319 std::unique_ptr<QgsPointCloudBlock> block( index->nodeData( n, request ) );
320 if ( !block )
321 continue;
322
323 const QgsVector3D blockScale = block->scale();
324 const QgsVector3D blockOffset = block->offset();
325
326 const char *ptr = block->data();
327 const QgsPointCloudAttributeCollection blockAttributes = block->attributes();
328 const std::size_t recordSize = blockAttributes.pointRecordSize();
329 int xOffset = 0, yOffset = 0, zOffset = 0;
330 const QgsPointCloudAttribute::DataType xType = blockAttributes.find( QStringLiteral( "X" ), xOffset )->type();
331 const QgsPointCloudAttribute::DataType yType = blockAttributes.find( QStringLiteral( "Y" ), yOffset )->type();
332 const QgsPointCloudAttribute::DataType zType = blockAttributes.find( QStringLiteral( "Z" ), zOffset )->type();
333 for ( int i = 0; i < block->pointCount(); ++i )
334 {
335 double x, y, z;
336 QgsPointCloudAttribute::getPointXYZ( ptr, i, recordSize, xOffset, xType, yOffset, yType, zOffset, zType, blockScale, blockOffset, x, y, z );
337 const QgsVector3D point( x, y, z );
338
339 // check whether point is in front of the ray
340 // similar to QgsRay3D::isInFront(), but using doubles
341 QgsVector3D vectorToPoint = point - adjustedRayOrigin;
342 vectorToPoint.normalize();
343 if ( QgsVector3D::dotProduct( vectorToPoint, adjustedRayDirection ) < 0.0 )
344 continue;
345
346 // calculate the angle between the point and the projected point
347 // similar to QgsRay3D::angleToPoint(), but using doubles
348 const QgsVector3D projPoint = adjustedRayOrigin + adjustedRayDirection * QgsVector3D::dotProduct( point - adjustedRayOrigin, adjustedRayDirection );
349 const QgsVector3D v1 = projPoint - adjustedRayOrigin ;
350 const QgsVector3D v2 = point - projPoint;
351 double angle = std::atan2( v2.length(), v1.length() ) * 180 / M_PI;
352 if ( angle > limitAngle )
353 continue;
354
355 const double dist = rayOriginMapCoords.distance( point );
356
357 if ( minDist < 0 || dist < minDist )
358 {
359 minDist = dist;
360 }
361 else if ( context.singleResult )
362 {
363 continue;
364 }
365
366 // Note : applying elevation properties is done in fromPointCloudIdentificationToIdentifyResults
367 QVariantMap pointAttr = QgsPointCloudAttribute::getAttributeMap( ptr, i * recordSize, blockAttributes );
368 pointAttr[ QStringLiteral( "X" ) ] = x;
369 pointAttr[ QStringLiteral( "Y" ) ] = y;
370 pointAttr[ QStringLiteral( "Z" ) ] = z;
371
372 const QgsVector3D worldPoint = factory->mRenderContext.mapToWorldCoordinates( point );
373 QgsRayCastingUtils::RayHit hit( dist, worldPoint.toVector3D(), FID_NULL, pointAttr );
374 if ( context.singleResult )
375 result.clear();
376 result.append( hit );
377 }
378 }
379 return result;
380}
Represents a indexed point cloud node in octree.
The Qgis class provides global constants for use throughout the application.
Definition qgis.h:54
@ Replacement
When tile is refined then its children should be used in place of itself.
@ Reverse
Reverse/inverse transform (from destination to source)
QgsVector3D origin() const
Returns coordinates in map CRS at which 3D scene has origin (0,0,0)
static QgsVector3D mapToWorldCoordinates(const QgsVector3D &mapCoords, const QgsVector3D &origin)
Converts map coordinates to 3D world coordinates (applies offset and turns (x,y,z) into (x,...
float yMax
Definition qgsaabb.h:90
float xMax
Definition qgsaabb.h:89
float xCenter() const
Returns center in X axis.
Definition qgsaabb.h:49
float xMin
Definition qgsaabb.h:86
float zMax
Definition qgsaabb.h:91
float yMin
Definition qgsaabb.h:87
float yCenter() const
Returns center in Y axis.
Definition qgsaabb.h:51
float zMin
Definition qgsaabb.h:88
float zCenter() const
Returns center in Z axis.
Definition qgsaabb.h:53
QgsPointCloudCategoryList getFilteredOutCategories() const
Gets the list of categories of the classification that should not be rendered.
Class for doing transforms between two map coordinate systems.
void setBallparkTransformsAreAppropriate(bool appropriate)
Sets whether approximate "ballpark" results are appropriate for this coordinate transform.
QgsPointXY transform(const QgsPointXY &point, Qgis::TransformDirection direction=Qgis::TransformDirection::Forward) const
Transform the point from the source CRS to the destination CRS.
Custom exception class for Coordinate Reference System related exceptions.
float pointSize() const
Returns the point size of the point cloud.
Collection of point cloud attributes.
int pointRecordSize() const
Returns total size of record.
const QgsPointCloudAttribute * find(const QString &attributeName, int &offset) const
Finds the attribute with the name.
QVector< QgsPointCloudAttribute > attributes() const
Returns all attributes.
DataType
Systems of unit measurement.
static void getPointXYZ(const char *ptr, int i, std::size_t pointRecordSize, int xOffset, QgsPointCloudAttribute::DataType xType, int yOffset, QgsPointCloudAttribute::DataType yType, int zOffset, QgsPointCloudAttribute::DataType zType, const QgsVector3D &indexScale, const QgsVector3D &indexOffset, double &x, double &y, double &z)
Retrieves the x, y, z values for the point at index i.
static QVariantMap getAttributeMap(const char *data, std::size_t recordOffset, const QgsPointCloudAttributeCollection &attributeCollection)
Retrieves all the attributes of a point.
DataType type() const
Returns the data type.
Represents packaged data bounds.
qint64 xMin() const
Returns x min.
qint64 zMin() const
Returns z min.
qint64 yMax() const
Returns y max.
qint64 xMax() const
Returns x max.
qint64 zMax() const
Returns z max.
qint64 yMin() const
Returns y min.
Represents a indexed point clouds data in octree.
virtual bool hasNode(const IndexedPointCloudNode &n) const
Returns whether the octree contain given node.
virtual std::unique_ptr< QgsPointCloudBlock > nodeData(const IndexedPointCloudNode &n, const QgsPointCloudRequest &request)=0
Returns node data block.
void setAttributes(const QgsPointCloudAttributeCollection &attributes)
Sets native attributes of the data.
QgsPointCloudAttributeCollection attributes() const
Returns all attributes that are stored in the file.
Point cloud data request.
void setAttributes(const QgsPointCloudAttributeCollection &attributes)
Set attributes filter in the request.
Class for storage of 3D vectors similar to QVector3D, with the difference that it uses double precisi...
Definition qgsvector3d.h:31
double y() const
Returns Y coordinate.
Definition qgsvector3d.h:50
double z() const
Returns Z coordinate.
Definition qgsvector3d.h:52
QVector3D toVector3D() const
Converts the current object to QVector3D.
static double dotProduct(const QgsVector3D &v1, const QgsVector3D &v2)
Returns the dot product of two vectors.
double x() const
Returns X coordinate.
Definition qgsvector3d.h:48
void normalize()
Normalizes the current vector in place.
double length() const
Returns the length of the vector.
double ANALYSIS_EXPORT angle(QgsPoint *p1, QgsPoint *p2, QgsPoint *p3, QgsPoint *p4)
Calculates the angle between two segments (in 2 dimension, z-values are ignored)
#define FID_NULL
#define QgsDebugMsgLevel(str, level)
Definition qgslogger.h:39
#define QgsDebugError(str)
Definition qgslogger.h:38
Helper struct to store ray casting parameters.
QSize screenSize
QSize of the 3d engine window.
bool singleResult
If set to true, only the closest point cloud hit will be returned (other entities always return only ...
Helper struct to store ray casting results.