QGIS API Documentation 3.29.0-Master (ade4f0cf0f)
qgseptpointcloudindex.cpp
Go to the documentation of this file.
1/***************************************************************************
2 qgspointcloudindex.cpp
3 --------------------
4 begin : October 2020
5 copyright : (C) 2020 by Peter Petrik
6 email : zilolv 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 ***************************************************************************/
17
19#include <QFile>
20#include <QFileInfo>
21#include <QDir>
22#include <QJsonArray>
23#include <QJsonDocument>
24#include <QJsonObject>
25#include <QTime>
26#include <QtDebug>
27#include <QQueue>
28
29#include "qgseptdecoder.h"
30#include "qgslazdecoder.h"
34#include "qgslogger.h"
35#include "qgsfeedback.h"
36#include "qgsmessagelog.h"
37#include "qgspointcloudexpression.h"
38
40
41#define PROVIDER_KEY QStringLiteral( "ept" )
42#define PROVIDER_DESCRIPTION QStringLiteral( "EPT point cloud provider" )
43
44QgsEptPointCloudIndex::QgsEptPointCloudIndex() = default;
45
46QgsEptPointCloudIndex::~QgsEptPointCloudIndex() = default;
47
48std::unique_ptr<QgsPointCloudIndex> QgsEptPointCloudIndex::clone() const
49{
50 QgsEptPointCloudIndex *clone = new QgsEptPointCloudIndex;
51 QMutexLocker locker( &mHierarchyMutex );
52 copyCommonProperties( clone );
53 return std::unique_ptr<QgsPointCloudIndex>( clone );
54}
55
56void QgsEptPointCloudIndex::load( const QString &fileName )
57{
58 QFile f( fileName );
59 if ( !f.open( QIODevice::ReadOnly ) )
60 {
61 mError = tr( "Unable to open %1 for reading" ).arg( fileName );
62 mIsValid = false;
63 return;
64 }
65
66 const QDir directory = QFileInfo( fileName ).absoluteDir();
67 mDirectory = directory.absolutePath();
68
69 const QByteArray dataJson = f.readAll();
70 bool success = loadSchema( dataJson );
71
72 if ( success )
73 {
74 // try to import the metadata too!
75 QFile manifestFile( mDirectory + QStringLiteral( "/ept-sources/manifest.json" ) );
76 if ( manifestFile.open( QIODevice::ReadOnly ) )
77 {
78 const QByteArray manifestJson = manifestFile.readAll();
79 loadManifest( manifestJson );
80 }
81 }
82
83 if ( success )
84 {
85 success = loadHierarchy();
86 }
87
88 mIsValid = success;
89}
90
91void QgsEptPointCloudIndex::loadManifest( const QByteArray &manifestJson )
92{
93 QJsonParseError err;
94 // try to import the metadata too!
95 const QJsonDocument manifestDoc = QJsonDocument::fromJson( manifestJson, &err );
96 if ( err.error == QJsonParseError::NoError )
97 {
98 const QJsonArray manifestArray = manifestDoc.array();
99 // TODO how to handle multiple?
100 if ( ! manifestArray.empty() )
101 {
102 const QJsonObject sourceObject = manifestArray.at( 0 ).toObject();
103 const QString metadataPath = sourceObject.value( QStringLiteral( "metadataPath" ) ).toString();
104 QFile metadataFile( mDirectory + QStringLiteral( "/ept-sources/" ) + metadataPath );
105 if ( metadataFile.open( QIODevice::ReadOnly ) )
106 {
107 const QByteArray metadataJson = metadataFile.readAll();
108 const QJsonDocument metadataDoc = QJsonDocument::fromJson( metadataJson, &err );
109 if ( err.error == QJsonParseError::NoError )
110 {
111 const QJsonObject metadataObject = metadataDoc.object().value( QStringLiteral( "metadata" ) ).toObject();
112 if ( !metadataObject.empty() )
113 {
114 const QJsonObject sourceMetadata = metadataObject.constBegin().value().toObject();
115 mOriginalMetadata = sourceMetadata.toVariantMap();
116 }
117 }
118 }
119 }
120 }
121}
122
123bool QgsEptPointCloudIndex::loadSchema( const QByteArray &dataJson )
124{
125 QJsonParseError err;
126 const QJsonDocument doc = QJsonDocument::fromJson( dataJson, &err );
127 if ( err.error != QJsonParseError::NoError )
128 return false;
129 const QJsonObject result = doc.object();
130 mDataType = result.value( QLatin1String( "dataType" ) ).toString(); // "binary" or "laszip"
131 if ( mDataType != QLatin1String( "laszip" ) && mDataType != QLatin1String( "binary" ) && mDataType != QLatin1String( "zstandard" ) )
132 return false;
133
134 const QString hierarchyType = result.value( QLatin1String( "hierarchyType" ) ).toString(); // "json" or "gzip"
135 if ( hierarchyType != QLatin1String( "json" ) )
136 return false;
137
138 mSpan = result.value( QLatin1String( "span" ) ).toInt();
139 mPointCount = result.value( QLatin1String( "points" ) ).toDouble();
140
141 // WKT
142 const QJsonObject srs = result.value( QLatin1String( "srs" ) ).toObject();
143 mWkt = srs.value( QLatin1String( "wkt" ) ).toString();
144
145 // rectangular
146 const QJsonArray bounds = result.value( QLatin1String( "bounds" ) ).toArray();
147 if ( bounds.size() != 6 )
148 return false;
149
150 const QJsonArray boundsConforming = result.value( QLatin1String( "boundsConforming" ) ).toArray();
151 if ( boundsConforming.size() != 6 )
152 return false;
153 mExtent.set( boundsConforming[0].toDouble(), boundsConforming[1].toDouble(),
154 boundsConforming[3].toDouble(), boundsConforming[4].toDouble() );
155 mZMin = boundsConforming[2].toDouble();
156 mZMax = boundsConforming[5].toDouble();
157
158 const QJsonArray schemaArray = result.value( QLatin1String( "schema" ) ).toArray();
160
161 for ( const QJsonValue &schemaItem : schemaArray )
162 {
163 const QJsonObject schemaObj = schemaItem.toObject();
164 const QString name = schemaObj.value( QLatin1String( "name" ) ).toString();
165 const QString type = schemaObj.value( QLatin1String( "type" ) ).toString();
166
167 const int size = schemaObj.value( QLatin1String( "size" ) ).toInt();
168
169 if ( type == QLatin1String( "float" ) && ( size == 4 ) )
170 {
172 }
173 else if ( type == QLatin1String( "float" ) && ( size == 8 ) )
174 {
176 }
177 else if ( size == 1 )
178 {
180 }
181 else if ( type == QLatin1String( "unsigned" ) && size == 2 )
182 {
184 }
185 else if ( size == 2 )
186 {
188 }
189 else if ( size == 4 )
190 {
192 }
193 else
194 {
195 // unknown attribute type
196 return false;
197 }
198
199 double scale = 1.f;
200 if ( schemaObj.contains( QLatin1String( "scale" ) ) )
201 scale = schemaObj.value( QLatin1String( "scale" ) ).toDouble();
202
203 double offset = 0.f;
204 if ( schemaObj.contains( QLatin1String( "offset" ) ) )
205 offset = schemaObj.value( QLatin1String( "offset" ) ).toDouble();
206
207 if ( name == QLatin1String( "X" ) )
208 {
209 mOffset.set( offset, mOffset.y(), mOffset.z() );
210 mScale.set( scale, mScale.y(), mScale.z() );
211 }
212 else if ( name == QLatin1String( "Y" ) )
213 {
214 mOffset.set( mOffset.x(), offset, mOffset.z() );
215 mScale.set( mScale.x(), scale, mScale.z() );
216 }
217 else if ( name == QLatin1String( "Z" ) )
218 {
219 mOffset.set( mOffset.x(), mOffset.y(), offset );
220 mScale.set( mScale.x(), mScale.y(), scale );
221 }
222
223 // store any metadata stats which are present for the attribute
224 AttributeStatistics stats;
225 bool foundStats = false;
226 if ( schemaObj.contains( QLatin1String( "count" ) ) )
227 {
228 stats.count = schemaObj.value( QLatin1String( "count" ) ).toInt();
229 foundStats = true;
230 }
231 if ( schemaObj.contains( QLatin1String( "minimum" ) ) )
232 {
233 stats.minimum = schemaObj.value( QLatin1String( "minimum" ) ).toDouble();
234 foundStats = true;
235 }
236 if ( schemaObj.contains( QLatin1String( "maximum" ) ) )
237 {
238 stats.maximum = schemaObj.value( QLatin1String( "maximum" ) ).toDouble();
239 foundStats = true;
240 }
241 if ( schemaObj.contains( QLatin1String( "count" ) ) )
242 {
243 stats.mean = schemaObj.value( QLatin1String( "mean" ) ).toDouble();
244 foundStats = true;
245 }
246 if ( schemaObj.contains( QLatin1String( "stddev" ) ) )
247 {
248 stats.stDev = schemaObj.value( QLatin1String( "stddev" ) ).toDouble();
249 foundStats = true;
250 }
251 if ( schemaObj.contains( QLatin1String( "variance" ) ) )
252 {
253 stats.variance = schemaObj.value( QLatin1String( "variance" ) ).toDouble();
254 foundStats = true;
255 }
256 if ( foundStats )
257 mMetadataStats.insert( name, stats );
258
259 if ( schemaObj.contains( QLatin1String( "counts" ) ) )
260 {
261 QMap< int, int > classCounts;
262 const QJsonArray counts = schemaObj.value( QLatin1String( "counts" ) ).toArray();
263 for ( const QJsonValue &count : counts )
264 {
265 const QJsonObject countObj = count.toObject();
266 classCounts.insert( countObj.value( QLatin1String( "value" ) ).toInt(), countObj.value( QLatin1String( "count" ) ).toInt() );
267 }
268 mAttributeClasses.insert( name, classCounts );
269 }
270 }
271 setAttributes( attributes );
272
273 // save mRootBounds
274
275 // bounds (cube - octree volume)
276 const double xmin = bounds[0].toDouble();
277 const double ymin = bounds[1].toDouble();
278 const double zmin = bounds[2].toDouble();
279 const double xmax = bounds[3].toDouble();
280 const double ymax = bounds[4].toDouble();
281 const double zmax = bounds[5].toDouble();
282
283 mRootBounds = QgsPointCloudDataBounds(
284 ( xmin - mOffset.x() ) / mScale.x(),
285 ( ymin - mOffset.y() ) / mScale.y(),
286 ( zmin - mOffset.z() ) / mScale.z(),
287 ( xmax - mOffset.x() ) / mScale.x(),
288 ( ymax - mOffset.y() ) / mScale.y(),
289 ( zmax - mOffset.z() ) / mScale.z()
290 );
291
292
293#ifdef QGIS_DEBUG
294 double dx = xmax - xmin, dy = ymax - ymin, dz = zmax - zmin;
295 QgsDebugMsgLevel( QStringLiteral( "lvl0 node size in CRS units: %1 %2 %3" ).arg( dx ).arg( dy ).arg( dz ), 2 ); // all dims should be the same
296 QgsDebugMsgLevel( QStringLiteral( "res at lvl0 %1" ).arg( dx / mSpan ), 2 );
297 QgsDebugMsgLevel( QStringLiteral( "res at lvl1 %1" ).arg( dx / mSpan / 2 ), 2 );
298 QgsDebugMsgLevel( QStringLiteral( "res at lvl2 %1 with node size %2" ).arg( dx / mSpan / 4 ).arg( dx / 4 ), 2 );
299#endif
300
301 return true;
302}
303
304QgsPointCloudBlock *QgsEptPointCloudIndex::nodeData( const IndexedPointCloudNode &n, const QgsPointCloudRequest &request )
305{
306 mHierarchyMutex.lock();
307 const bool found = mHierarchy.contains( n );
308 mHierarchyMutex.unlock();
309 if ( !found )
310 return nullptr;
311
312 // we need to create a copy of the expression to pass to the decoder
313 // as the same QgsPointCloudExpression object mighgt be concurrently
314 // used on another thread, for example in a 3d view
315 QgsPointCloudExpression filterExpression = mFilterExpression;
316 QgsPointCloudAttributeCollection requestAttributes = request.attributes();
317 requestAttributes.extend( attributes(), filterExpression.referencedAttributes() );
318
319 if ( mDataType == QLatin1String( "binary" ) )
320 {
321 const QString filename = QStringLiteral( "%1/ept-data/%2.bin" ).arg( mDirectory, n.toString() );
322 return QgsEptDecoder::decompressBinary( filename, attributes(), requestAttributes, scale(), offset(), filterExpression );
323 }
324 else if ( mDataType == QLatin1String( "zstandard" ) )
325 {
326 const QString filename = QStringLiteral( "%1/ept-data/%2.zst" ).arg( mDirectory, n.toString() );
327 return QgsEptDecoder::decompressZStandard( filename, attributes(), request.attributes(), scale(), offset(), filterExpression );
328 }
329 else if ( mDataType == QLatin1String( "laszip" ) )
330 {
331 const QString filename = QStringLiteral( "%1/ept-data/%2.laz" ).arg( mDirectory, n.toString() );
332 return QgsLazDecoder::decompressLaz( filename, requestAttributes, filterExpression );
333 }
334 else
335 {
336 return nullptr; // unsupported
337 }
338}
339
340QgsPointCloudBlockRequest *QgsEptPointCloudIndex::asyncNodeData( const IndexedPointCloudNode &n, const QgsPointCloudRequest &request )
341{
342 Q_UNUSED( n );
343 Q_UNUSED( request );
344 Q_ASSERT( false );
345 return nullptr; // unsupported
346}
347
349{
351}
352
353qint64 QgsEptPointCloudIndex::pointCount() const
354{
355 return mPointCount;
356}
357
358bool QgsEptPointCloudIndex::hasStatisticsMetadata() const
359{
360 return !mMetadataStats.isEmpty();
361}
362
363QVariant QgsEptPointCloudIndex::metadataStatistic( const QString &attribute, QgsStatisticalSummary::Statistic statistic ) const
364{
365 if ( !mMetadataStats.contains( attribute ) )
366 return QVariant();
367
368 const AttributeStatistics &stats = mMetadataStats[ attribute ];
369 switch ( statistic )
370 {
372 return stats.count >= 0 ? QVariant( stats.count ) : QVariant();
373
375 return std::isnan( stats.mean ) ? QVariant() : QVariant( stats.mean );
376
378 return std::isnan( stats.stDev ) ? QVariant() : QVariant( stats.stDev );
379
381 return stats.minimum;
382
384 return stats.maximum;
385
387 return stats.minimum.isValid() && stats.maximum.isValid() ? QVariant( stats.maximum.toDouble() - stats.minimum.toDouble() ) : QVariant();
388
402 return QVariant();
403 }
404 return QVariant();
405}
406
407QVariantList QgsEptPointCloudIndex::metadataClasses( const QString &attribute ) const
408{
409 QVariantList classes;
410 const QMap< int, int > values = mAttributeClasses.value( attribute );
411 for ( auto it = values.constBegin(); it != values.constEnd(); ++it )
412 {
413 classes << it.key();
414 }
415 return classes;
416}
417
418QVariant QgsEptPointCloudIndex::metadataClassStatistic( const QString &attribute, const QVariant &value, QgsStatisticalSummary::Statistic statistic ) const
419{
420 if ( statistic != QgsStatisticalSummary::Count )
421 return QVariant();
422
423 const QMap< int, int > values = mAttributeClasses.value( attribute );
424 if ( !values.contains( value.toInt() ) )
425 return QVariant();
426 return values.value( value.toInt() );
427}
428
429bool QgsEptPointCloudIndex::loadHierarchy()
430{
431 QQueue<QString> queue;
432 queue.enqueue( QStringLiteral( "0-0-0-0" ) );
433 while ( !queue.isEmpty() )
434 {
435 const QString filename = QStringLiteral( "%1/ept-hierarchy/%2.json" ).arg( mDirectory, queue.dequeue() );
436 QFile fH( filename );
437 if ( !fH.open( QIODevice::ReadOnly ) )
438 {
439 QgsDebugMsgLevel( QStringLiteral( "unable to read hierarchy from file %1" ).arg( filename ), 2 );
440 mError = QStringLiteral( "unable to read hierarchy from file %1" ).arg( filename );
441 return false;
442 }
443
444 const QByteArray dataJsonH = fH.readAll();
445 QJsonParseError errH;
446 const QJsonDocument docH = QJsonDocument::fromJson( dataJsonH, &errH );
447 if ( errH.error != QJsonParseError::NoError )
448 {
449 QgsDebugMsgLevel( QStringLiteral( "QJsonParseError when reading hierarchy from file %1" ).arg( filename ), 2 );
450 mError = QStringLiteral( "QJsonParseError when reading hierarchy from file %1" ).arg( filename );
451 return false;
452 }
453
454 const QJsonObject rootHObj = docH.object();
455 for ( auto it = rootHObj.constBegin(); it != rootHObj.constEnd(); ++it )
456 {
457 const QString nodeIdStr = it.key();
458 const int nodePointCount = it.value().toInt();
459 if ( nodePointCount < 0 )
460 {
461 queue.enqueue( nodeIdStr );
462 }
463 else
464 {
466 mHierarchyMutex.lock();
467 mHierarchy[nodeId] = nodePointCount;
468 mHierarchyMutex.unlock();
469 }
470 }
471 }
472 return true;
473}
474
475bool QgsEptPointCloudIndex::isValid() const
476{
477 return mIsValid;
478}
479
480void QgsEptPointCloudIndex::copyCommonProperties( QgsEptPointCloudIndex *destination ) const
481{
483
484 // QgsEptPointCloudIndex specific fields
485 destination->mIsValid = mIsValid;
486 destination->mDataType = mDataType;
487 destination->mDirectory = mDirectory;
488 destination->mWkt = mWkt;
489 destination->mPointCount = mPointCount;
490 destination->mMetadataStats = mMetadataStats;
491 destination->mAttributeClasses = mAttributeClasses;
492 destination->mOriginalMetadata = mOriginalMetadata;
493}
494
Represents a indexed point cloud node in octree.
static IndexedPointCloudNode fromString(const QString &str)
Creates node from string.
QString toString() const
Encode node to string.
This class represents a coordinate reference system (CRS).
static QgsCoordinateReferenceSystem fromWkt(const QString &wkt)
Creates a CRS from a WKT spatial ref sys definition string.
Collection of point cloud attributes.
void push_back(const QgsPointCloudAttribute &attribute)
Adds extra attribute.
void extend(const QgsPointCloudAttributeCollection &otherCollection, const QSet< QString > &matchingNames)
Adds specific missing attributes from another QgsPointCloudAttributeCollection.
Attribute for point cloud data pair of name and size in bytes.
@ UShort
Unsigned short int 2 bytes.
Base class for handling loading QgsPointCloudBlock asynchronously.
Base class for storing raw data from point cloud nodes.
Represents packaged data bounds.
void copyCommonProperties(QgsPointCloudIndex *destination) const
Copies common properties to the destination index.
Point cloud data request.
QgsPointCloudAttributeCollection attributes() const
Returns attributes.
Statistic
Enumeration of flags that specify statistics to be calculated.
@ InterQuartileRange
Inter quartile range (IQR)
@ Median
Median of values.
@ Variety
Variety (count of distinct) values.
@ First
First value (since QGIS 3.6)
@ Last
Last value (since QGIS 3.6)
@ Minority
Minority of values.
@ ThirdQuartile
Third quartile.
@ Majority
Majority of values.
@ CountMissing
Number of missing (null) values.
@ Range
Range of values (max - min)
@ StDevSample
Sample standard deviation of values.
@ FirstQuartile
First quartile.
@ StDev
Standard deviation of values.
#define QgsDebugMsgLevel(str, level)
Definition: qgslogger.h:39
const QgsCoordinateReferenceSystem & crs