QGIS API Documentation 3.38.0-Grenoble (exported)
Loading...
Searching...
No Matches
qgsraycastingutils.cpp
Go to the documentation of this file.
1/***************************************************************************
2 qgsraycastingutils_h.cpp
3 --------------------------------------
4 Date : June 2018
5 Copyright : (C) 2018 by Martin Dobias
6 Email : wonder 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 "qgis.h"
19#include "qgsaabb.h"
20#include "qgslogger.h"
21
22#include <Qt3DRender/QCamera>
23#include <Qt3DRender/QGeometryRenderer>
24
25#if QT_VERSION < QT_VERSION_CHECK(6, 0, 0)
26#include <Qt3DRender/QAttribute>
27#include <Qt3DRender/QBuffer>
28#include <Qt3DRender/QGeometry>
29typedef Qt3DRender::QAttribute Qt3DQAttribute;
30typedef Qt3DRender::QBuffer Qt3DQBuffer;
31typedef Qt3DRender::QGeometry Qt3DQGeometry;
32#else
33#include <Qt3DCore/QAttribute>
34#include <Qt3DCore/QBuffer>
35#include <Qt3DCore/QGeometry>
36typedef Qt3DCore::QAttribute Qt3DQAttribute;
37typedef Qt3DCore::QBuffer Qt3DQBuffer;
38typedef Qt3DCore::QGeometry Qt3DQGeometry;
39#endif
40
42
43
44namespace QgsRayCastingUtils
45{
46
47// copied from qt3d/src/render/raycasting/qray3d_p.h + qray3d.cpp
48// by KDAB, licensed under the terms of LGPL
49
50
51 Ray3D::Ray3D()
52 : m_direction( 0.0f, 0.0f, 1.0f )
53 {
54 }
55
56 Ray3D::Ray3D( QVector3D origin, QVector3D direction, float distance )
57 : m_origin( origin )
58 , m_direction( direction )
59 , m_distance( distance )
60 {}
61
62 QVector3D Ray3D::origin() const
63 {
64 return m_origin;
65 }
66
67 void Ray3D::setOrigin( QVector3D value )
68 {
69 m_origin = value;
70 }
71
72 QVector3D Ray3D::direction() const
73 {
74 return m_direction;
75 }
76
77 void Ray3D::setDirection( QVector3D value )
78 {
79 if ( value.isNull() )
80 return;
81
82 m_direction = value;
83 }
84
85 float Ray3D::distance() const
86 {
87 return m_distance;
88 }
89
90 void Ray3D::setDistance( float distance )
91 {
92 m_distance = distance;
93 }
94
95 QVector3D Ray3D::point( float t ) const
96 {
97 return m_origin + t * m_direction;
98 }
99
100 Ray3D &Ray3D::transform( const QMatrix4x4 &matrix )
101 {
102 m_origin = matrix * m_origin;
103 m_direction = matrix.mapVector( m_direction );
104
105 return *this;
106 }
107
108 Ray3D Ray3D::transformed( const QMatrix4x4 &matrix ) const
109 {
110 return Ray3D( matrix * m_origin, matrix.mapVector( m_direction ) );
111 }
112
113 bool Ray3D::operator==( const Ray3D &other ) const
114 {
115 return m_origin == other.origin() && m_direction == other.direction();
116 }
117
118 bool Ray3D::operator!=( const Ray3D &other ) const
119 {
120 return !( *this == other );
121 }
122
123 bool Ray3D::contains( QVector3D point ) const
124 {
125 const QVector3D ppVec( point - m_origin );
126 if ( ppVec.isNull() ) // point coincides with origin
127 return true;
128 const float dot = QVector3D::dotProduct( ppVec, m_direction );
129 if ( qFuzzyIsNull( dot ) )
130 return false;
131 return qFuzzyCompare( dot * dot, ppVec.lengthSquared() * m_direction.lengthSquared() );
132 }
133
134 bool Ray3D::contains( const Ray3D &ray ) const
135 {
136 const float dot = QVector3D::dotProduct( m_direction, ray.direction() );
137 if ( !qFuzzyCompare( dot * dot, m_direction.lengthSquared() * ray.direction().lengthSquared() ) )
138 return false;
139 return contains( ray.origin() );
140 }
141
142 float Ray3D::projectedDistance( QVector3D point ) const
143 {
144 Q_ASSERT( !m_direction.isNull() );
145
146 return QVector3D::dotProduct( point - m_origin, m_direction ) /
147 m_direction.lengthSquared();
148 }
149
150 QVector3D Ray3D::project( QVector3D vector ) const
151 {
152 const QVector3D norm = m_direction.normalized();
153 return QVector3D::dotProduct( vector, norm ) * norm;
154 }
155
156
157 float Ray3D::distance( QVector3D point ) const
158 {
159 const float t = projectedDistance( point );
160 return ( point - ( m_origin + t * m_direction ) ).length();
161 }
162
163 QDebug operator<<( QDebug dbg, const Ray3D &ray )
164 {
165 const QDebugStateSaver saver( dbg );
166 dbg.nospace() << "QRay3D(origin("
167 << ray.origin().x() << ", " << ray.origin().y() << ", "
168 << ray.origin().z() << ") - direction("
169 << ray.direction().x() << ", " << ray.direction().y() << ", "
170 << ray.direction().z() << "))";
171 return dbg;
172 }
173
174}
175
176
178
179
180struct box
181{
182 box( const QgsAABB &b )
183 {
184 min[0] = b.xMin; min[1] = b.yMin; min[2] = b.zMin;
185 max[0] = b.xMax; max[1] = b.yMax; max[2] = b.zMax;
186 }
187 double min[3];
188 double max[3];
189};
190
191struct ray
192{
193 ray( double xO, double yO, double zO, double xD, double yD, double zD )
194 {
195 origin[0] = xO; origin[1] = yO; origin[2] = zO;
196 dir[0] = xD; dir[1] = yD; dir[2] = zD;
197 dir_inv[0] = 1 / dir[0]; dir_inv[1] = 1 / dir[1]; dir_inv[2] = 1 / dir[2];
198 }
199 ray( const QgsRayCastingUtils::Ray3D &r )
200 {
201 // ignoring length...
202 origin[0] = r.origin().x(); origin[1] = r.origin().y(); origin[2] = r.origin().z();
203 dir[0] = r.direction().x(); dir[1] = r.direction().y(); dir[2] = r.direction().z();
204 dir_inv[0] = 1 / dir[0]; dir_inv[1] = 1 / dir[1]; dir_inv[2] = 1 / dir[2];
205 }
206 double origin[3];
207 double dir[3];
208 double dir_inv[3];
209};
210
211// https://tavianator.com/fast-branchless-raybounding-box-intersections/
212// https://tavianator.com/fast-branchless-raybounding-box-intersections-part-2-nans/
213
214bool intersection( const box &b, const ray &r )
215{
216 double t1 = ( b.min[0] - r.origin[0] ) * r.dir_inv[0];
217 double t2 = ( b.max[0] - r.origin[0] ) * r.dir_inv[0];
218
219 double tmin = std::min( t1, t2 );
220 double tmax = std::max( t1, t2 );
221
222 for ( int i = 1; i < 3; ++i )
223 {
224 t1 = ( b.min[i] - r.origin[i] ) * r.dir_inv[i];
225 t2 = ( b.max[i] - r.origin[i] ) * r.dir_inv[i];
226
227 tmin = std::max( tmin, std::min( std::min( t1, t2 ), tmax ) );
228 tmax = std::min( tmax, std::max( std::max( t1, t2 ), tmin ) );
229 }
230
231 return tmax > std::max( tmin, 0.0 );
232}
233
234
235namespace QgsRayCastingUtils
236{
237
238 bool rayBoxIntersection( const Ray3D &r, const QgsAABB &aabb )
239 {
240 box b( aabb );
241
242 // intersection() does not like yMin==yMax (excludes borders)
243 if ( b.min[0] == b.max[0] ) b.max[0] += 0.1;
244 if ( b.min[1] == b.max[1] ) b.max[1] += 0.1;
245 if ( b.min[2] == b.max[2] ) b.max[2] += 0.1;
246
247 return intersection( b, ray( r ) );
248 }
249
250// copied from intersectsSegmentTriangle() from qt3d/src/render/backend/triangleboundingvolume.cpp
251// by KDAB, licensed under the terms of LGPL
252 bool rayTriangleIntersection( const Ray3D &ray,
253 QVector3D a,
254 QVector3D b,
255 QVector3D c,
256 QVector3D &uvw,
257 float &t )
258 {
259 // Note: a, b, c in clockwise order
260 // RealTime Collision Detection page 192
261
262 const QVector3D ab = b - a;
263 const QVector3D ac = c - a;
264 const QVector3D qp = ( ray.origin() - ray.point( ray.distance() ) );
265
266 const QVector3D n = QVector3D::crossProduct( ab, ac );
267 const float d = QVector3D::dotProduct( qp, n );
268
269 if ( d <= 0.0f || std::isnan( d ) )
270 return false;
271
272 const QVector3D ap = ray.origin() - a;
273 t = QVector3D::dotProduct( ap, n );
274
275 if ( t < 0.0f || t > d )
276 return false;
277
278 const QVector3D e = QVector3D::crossProduct( qp, ap );
279 uvw.setY( QVector3D::dotProduct( ac, e ) );
280
281 if ( uvw.y() < 0.0f || uvw.y() > d )
282 return false;
283
284 uvw.setZ( -QVector3D::dotProduct( ab, e ) );
285
286 if ( uvw.z() < 0.0f || uvw.y() + uvw.z() > d )
287 return false;
288
289 const float ood = 1.0f / d;
290 t *= ood;
291 uvw.setY( uvw.y() * ood );
292 uvw.setZ( uvw.z() * ood );
293 uvw.setX( 1.0f - uvw.y() - uvw.z() );
294
295 return true;
296 }
297
298 bool rayMeshIntersection( Qt3DRender::QGeometryRenderer *geometryRenderer,
299 const QgsRayCastingUtils::Ray3D &r,
300 const QMatrix4x4 &worldTransform,
301 QVector3D &intPt,
302 int &triangleIndex )
303 {
304 if ( geometryRenderer->primitiveType() != Qt3DRender::QGeometryRenderer::Triangles )
305 {
306 QgsDebugError( QString( "Unsupported primitive type for intersection: " ).arg( geometryRenderer->primitiveType() ) );
307 return false;
308 }
309 if ( geometryRenderer->instanceCount() != 1 || geometryRenderer->indexOffset() != 0 || geometryRenderer->indexBufferByteOffset() != 0 || geometryRenderer->firstVertex() != 0 || geometryRenderer->firstInstance() != 0 )
310 {
311 QgsDebugError( QString( "Unsupported geometry renderer for intersection." ) );
312 return false;
313 }
314
315 Qt3DQGeometry *geometry = geometryRenderer->geometry();
316
317 Qt3DQAttribute *positionAttr = nullptr;
318 Qt3DQAttribute *indexAttr = nullptr;
319 for ( Qt3DQAttribute *attr : geometry->attributes() )
320 {
321 if ( attr->name() == Qt3DQAttribute::defaultPositionAttributeName() )
322 {
323 positionAttr = attr;
324 }
325 else if ( attr->attributeType() == Qt3DQAttribute::IndexAttribute )
326 {
327 indexAttr = attr;
328 }
329 }
330
331 if ( !positionAttr )
332 {
333 QgsDebugError( "Could not find position attribute!" );
334 return false;
335 }
336
337 if ( positionAttr->vertexBaseType() != Qt3DQAttribute::Float || positionAttr->vertexSize() != 3 )
338 {
339 QgsDebugError( QString( "Unsupported position attribute: base type %1, vertex size %2" ). arg( positionAttr->vertexBaseType() ).arg( positionAttr->vertexSize() ) );
340 return false;
341 }
342
343 const QByteArray vertexBuf = positionAttr->buffer()->data();
344 const char *vertexPtr = vertexBuf.constData();
345 vertexPtr += positionAttr->byteOffset();
346 int vertexByteStride = positionAttr->byteStride() == 0 ? 3 * sizeof( float ) : positionAttr->byteStride();
347
348 const uchar *indexPtrUChar = nullptr;
349 const ushort *indexPtrUShort = nullptr;
350 const uint *indexPtrUInt = nullptr;
351 if ( indexAttr )
352 {
353 if ( indexAttr->byteStride() != 0 || indexAttr->vertexSize() != 1 )
354 {
355 QgsDebugError( QString( "Unsupported index attribute: stride %1, vertex size %2" ).arg( indexAttr->byteStride() ).arg( indexAttr->vertexSize() ) );
356 return false;
357 }
358
359 const QByteArray indexBuf = indexAttr->buffer()->data();
360 if ( indexAttr->vertexBaseType() == Qt3DQAttribute::UnsignedByte )
361 {
362 indexPtrUChar = reinterpret_cast<const uchar *>( indexBuf.constData() + indexAttr->byteOffset() );
363 }
364 else if ( indexAttr->vertexBaseType() == Qt3DQAttribute::UnsignedShort )
365 {
366 indexPtrUShort = reinterpret_cast<const ushort *>( indexBuf.constData() + indexAttr->byteOffset() );
367 }
368 else if ( indexAttr->vertexBaseType() == Qt3DQAttribute::UnsignedInt )
369 {
370 indexPtrUInt = reinterpret_cast<const uint *>( indexBuf.constData() + indexAttr->byteOffset() );
371 }
372 else
373 {
374 QgsDebugError( QString( "Unsupported index attribute: base type %1" ).arg( indexAttr->vertexBaseType() ) );
375 return false;
376 }
377 }
378
379 int vertexCount = geometryRenderer->vertexCount();
380 if ( vertexCount == 0 && indexAttr )
381 {
382 vertexCount = indexAttr->count();
383 }
384 if ( vertexCount == 0 )
385 {
386 vertexCount = positionAttr->count();
387 }
388
389 QVector3D intersectionPt, minIntersectionPt;
390 float minDistance = -1;
391
392 for ( int i = 0; i < vertexCount; i += 3 )
393 {
394 int v0index = 0, v1index = 0, v2index = 0;
395 if ( !indexAttr )
396 {
397 v0index = i;
398 v1index = i + 1;
399 v2index = i + 2;
400 }
401 else if ( indexPtrUShort )
402 {
403 v0index = indexPtrUShort[i];
404 v1index = indexPtrUShort[i + 1];
405 v2index = indexPtrUShort[i + 2];
406 }
407 else if ( indexPtrUChar )
408 {
409 v0index = indexPtrUChar[i];
410 v1index = indexPtrUChar[i + 1];
411 v2index = indexPtrUChar[i + 2];
412 }
413 else if ( indexPtrUInt )
414 {
415 v0index = indexPtrUInt[i];
416 v1index = indexPtrUInt[i + 1];
417 v2index = indexPtrUInt[i + 2];
418 }
419 else
420 Q_ASSERT( false );
421
422 const float *v0ptr = reinterpret_cast<const float *>( vertexPtr + v0index * vertexByteStride );
423 const float *v1ptr = reinterpret_cast<const float *>( vertexPtr + v1index * vertexByteStride );
424 const float *v2ptr = reinterpret_cast<const float *>( vertexPtr + v2index * vertexByteStride );
425
426 const QVector3D a( v0ptr[0], v0ptr[1], v0ptr[2] );
427 const QVector3D b( v1ptr[0], v1ptr[1], v1ptr[2] );
428 const QVector3D c( v2ptr[0], v2ptr[1], v2ptr[2] );
429
430 // Currently the worldTransform only has vertical offset, so this could be optimized by applying the transform
431 // to the ray and the resulting intersecting point instead of all triangles
432 // Need to check for potential performance gains.
433 const QVector3D tA = worldTransform * a;
434 const QVector3D tB = worldTransform * b;
435 const QVector3D tC = worldTransform * c;
436
437 QVector3D uvw;
438 float t = 0;
439
440 // We're testing both triangle orientations here and ignoring the culling mode.
441 // We should probably respect the culling mode used for the entity and perform a
442 // single test using the properly oriented triangle.
443 if ( QgsRayCastingUtils::rayTriangleIntersection( r, tA, tB, tC, uvw, t ) ||
444 QgsRayCastingUtils::rayTriangleIntersection( r, tA, tC, tB, uvw, t ) )
445 {
446 intersectionPt = r.point( t * r.distance() );
447 const float distance = r.projectedDistance( intersectionPt );
448
449 // we only want the first intersection of the ray with the mesh (closest to the ray origin)
450 if ( minDistance == -1 || distance < minDistance )
451 {
452 triangleIndex = static_cast<int>( i / 3 );
453 minDistance = distance;
454 minIntersectionPt = intersectionPt;
455 }
456 }
457 }
458
459 if ( minDistance != -1 )
460 {
461 intPt = minIntersectionPt;
462 return true;
463 }
464 else
465 return false;
466 }
467}
468
469
float yMax
Definition qgsaabb.h:90
float xMax
Definition qgsaabb.h:89
float xMin
Definition qgsaabb.h:86
float zMax
Definition qgsaabb.h:91
float yMin
Definition qgsaabb.h:87
float zMin
Definition qgsaabb.h:88
As part of the API refactoring and improvements which landed in the Processing API was substantially reworked from the x version This was done in order to allow much of the underlying Processing framework to be ported into c
Qt3DCore::QAttribute Qt3DQAttribute
Qt3DCore::QGeometry Qt3DQGeometry
QDataStream & operator<<(QDataStream &out, const QgsFeature &feature)
Writes the feature to stream out. QGIS version compatibility is not guaranteed.
#define QgsDebugError(str)
Definition qgslogger.h:38
Qt3DCore::QAttribute Qt3DQAttribute
Qt3DCore::QBuffer Qt3DQBuffer
Qt3DCore::QGeometry Qt3DQGeometry