Branch data Line data Source code
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 : :
18 : : #include "qgseptpointcloudindex.h"
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 "qgscoordinatereferencesystem.h"
31 : : #include "qgspointcloudrequest.h"
32 : : #include "qgspointcloudattribute.h"
33 : : #include "qgslogger.h"
34 : : #include "qgsfeedback.h"
35 : : #include "qgsmessagelog.h"
36 : :
37 : : ///@cond PRIVATE
38 : :
39 : : #define PROVIDER_KEY QStringLiteral( "ept" )
40 : : #define PROVIDER_DESCRIPTION QStringLiteral( "EPT point cloud provider" )
41 : :
42 : 0 : QgsEptPointCloudIndex::QgsEptPointCloudIndex() = default;
43 : :
44 : 0 : QgsEptPointCloudIndex::~QgsEptPointCloudIndex() = default;
45 : :
46 : 0 : void QgsEptPointCloudIndex::load( const QString &fileName )
47 : : {
48 : 0 : QFile f( fileName );
49 : 0 : if ( !f.open( QIODevice::ReadOnly ) )
50 : : {
51 : 0 : QgsMessageLog::logMessage( tr( "Unable to open %1 for reading" ).arg( fileName ) );
52 : 0 : mIsValid = false;
53 : 0 : return;
54 : : }
55 : :
56 : 0 : const QDir directory = QFileInfo( fileName ).absoluteDir();
57 : 0 : mDirectory = directory.absolutePath();
58 : :
59 : 0 : const QByteArray dataJson = f.readAll();
60 : 0 : bool success = loadSchema( dataJson );
61 : :
62 : 0 : if ( success )
63 : : {
64 : : // try to import the metadata too!
65 : 0 : QFile manifestFile( mDirectory + QStringLiteral( "/ept-sources/manifest.json" ) );
66 : 0 : if ( manifestFile.open( QIODevice::ReadOnly ) )
67 : : {
68 : 0 : const QByteArray manifestJson = manifestFile.readAll();
69 : 0 : loadManifest( manifestJson );
70 : 0 : }
71 : 0 : }
72 : :
73 : 0 : if ( success )
74 : : {
75 : 0 : success = loadHierarchy();
76 : 0 : }
77 : :
78 : 0 : mIsValid = success;
79 : 0 : }
80 : :
81 : 0 : void QgsEptPointCloudIndex::loadManifest( const QByteArray &manifestJson )
82 : : {
83 : : QJsonParseError err;
84 : : // try to import the metadata too!
85 : 0 : const QJsonDocument manifestDoc = QJsonDocument::fromJson( manifestJson, &err );
86 : 0 : if ( err.error == QJsonParseError::NoError )
87 : : {
88 : 0 : const QJsonArray manifestArray = manifestDoc.array();
89 : : // TODO how to handle multiple?
90 : 0 : if ( ! manifestArray.empty() )
91 : : {
92 : 0 : const QJsonObject sourceObject = manifestArray.at( 0 ).toObject();
93 : 0 : const QString metadataPath = sourceObject.value( QStringLiteral( "metadataPath" ) ).toString();
94 : 0 : QFile metadataFile( mDirectory + QStringLiteral( "/ept-sources/" ) + metadataPath );
95 : 0 : if ( metadataFile.open( QIODevice::ReadOnly ) )
96 : : {
97 : 0 : const QByteArray metadataJson = metadataFile.readAll();
98 : 0 : const QJsonDocument metadataDoc = QJsonDocument::fromJson( metadataJson, &err );
99 : 0 : if ( err.error == QJsonParseError::NoError )
100 : : {
101 : 0 : const QJsonObject metadataObject = metadataDoc.object().value( QStringLiteral( "metadata" ) ).toObject();
102 : 0 : if ( !metadataObject.empty() )
103 : : {
104 : 0 : const QJsonObject sourceMetadata = metadataObject.constBegin().value().toObject();
105 : 0 : mOriginalMetadata = sourceMetadata.toVariantMap();
106 : 0 : }
107 : 0 : }
108 : 0 : }
109 : 0 : }
110 : 0 : }
111 : 0 : }
112 : :
113 : 0 : bool QgsEptPointCloudIndex::loadSchema( const QByteArray &dataJson )
114 : : {
115 : : QJsonParseError err;
116 : 0 : const QJsonDocument doc = QJsonDocument::fromJson( dataJson, &err );
117 : 0 : if ( err.error != QJsonParseError::NoError )
118 : 0 : return false;
119 : 0 : const QJsonObject result = doc.object();
120 : 0 : mDataType = result.value( QLatin1String( "dataType" ) ).toString(); // "binary" or "laszip"
121 : 0 : if ( mDataType != QLatin1String( "laszip" ) && mDataType != QLatin1String( "binary" ) && mDataType != QLatin1String( "zstandard" ) )
122 : 0 : return false;
123 : :
124 : 0 : const QString hierarchyType = result.value( QLatin1String( "hierarchyType" ) ).toString(); // "json" or "gzip"
125 : 0 : if ( hierarchyType != QLatin1String( "json" ) )
126 : 0 : return false;
127 : :
128 : 0 : mSpan = result.value( QLatin1String( "span" ) ).toInt();
129 : 0 : mPointCount = result.value( QLatin1String( "points" ) ).toDouble();
130 : :
131 : : // WKT
132 : 0 : const QJsonObject srs = result.value( QLatin1String( "srs" ) ).toObject();
133 : 0 : mWkt = srs.value( QLatin1String( "wkt" ) ).toString();
134 : :
135 : : // rectangular
136 : 0 : const QJsonArray bounds = result.value( QLatin1String( "bounds" ) ).toArray();
137 : 0 : if ( bounds.size() != 6 )
138 : 0 : return false;
139 : :
140 : 0 : const QJsonArray boundsConforming = result.value( QLatin1String( "boundsConforming" ) ).toArray();
141 : 0 : if ( boundsConforming.size() != 6 )
142 : 0 : return false;
143 : 0 : mExtent.set( boundsConforming[0].toDouble(), boundsConforming[1].toDouble(),
144 : 0 : boundsConforming[3].toDouble(), boundsConforming[4].toDouble() );
145 : 0 : mZMin = boundsConforming[2].toDouble();
146 : 0 : mZMax = boundsConforming[5].toDouble();
147 : :
148 : 0 : const QJsonArray schemaArray = result.value( QLatin1String( "schema" ) ).toArray();
149 : 0 : QgsPointCloudAttributeCollection attributes;
150 : :
151 : 0 : for ( const QJsonValue &schemaItem : schemaArray )
152 : : {
153 : 0 : const QJsonObject schemaObj = schemaItem.toObject();
154 : 0 : const QString name = schemaObj.value( QLatin1String( "name" ) ).toString();
155 : 0 : const QString type = schemaObj.value( QLatin1String( "type" ) ).toString();
156 : :
157 : 0 : int size = schemaObj.value( QLatin1String( "size" ) ).toInt();
158 : :
159 : 0 : if ( type == QLatin1String( "float" ) && ( size == 4 ) )
160 : : {
161 : 0 : attributes.push_back( QgsPointCloudAttribute( name, QgsPointCloudAttribute::Float ) );
162 : 0 : }
163 : 0 : else if ( type == QLatin1String( "float" ) && ( size == 8 ) )
164 : : {
165 : 0 : attributes.push_back( QgsPointCloudAttribute( name, QgsPointCloudAttribute::Double ) );
166 : 0 : }
167 : 0 : else if ( size == 1 )
168 : : {
169 : 0 : attributes.push_back( QgsPointCloudAttribute( name, QgsPointCloudAttribute::Char ) );
170 : 0 : }
171 : 0 : else if ( type == QLatin1String( "unsigned" ) && size == 2 )
172 : : {
173 : 0 : attributes.push_back( QgsPointCloudAttribute( name, QgsPointCloudAttribute::UShort ) );
174 : 0 : }
175 : 0 : else if ( size == 2 )
176 : : {
177 : 0 : attributes.push_back( QgsPointCloudAttribute( name, QgsPointCloudAttribute::Short ) );
178 : 0 : }
179 : 0 : else if ( size == 4 )
180 : : {
181 : 0 : attributes.push_back( QgsPointCloudAttribute( name, QgsPointCloudAttribute::Int32 ) );
182 : 0 : }
183 : : else
184 : : {
185 : : // unknown attribute type
186 : 0 : return false;
187 : : }
188 : :
189 : 0 : double scale = 1.f;
190 : 0 : if ( schemaObj.contains( QLatin1String( "scale" ) ) )
191 : 0 : scale = schemaObj.value( QLatin1String( "scale" ) ).toDouble();
192 : :
193 : 0 : double offset = 0.f;
194 : 0 : if ( schemaObj.contains( QLatin1String( "offset" ) ) )
195 : 0 : offset = schemaObj.value( QLatin1String( "offset" ) ).toDouble();
196 : :
197 : 0 : if ( name == QLatin1String( "X" ) )
198 : : {
199 : 0 : mOffset.set( offset, mOffset.y(), mOffset.z() );
200 : 0 : mScale.set( scale, mScale.y(), mScale.z() );
201 : 0 : }
202 : 0 : else if ( name == QLatin1String( "Y" ) )
203 : : {
204 : 0 : mOffset.set( mOffset.x(), offset, mOffset.z() );
205 : 0 : mScale.set( mScale.x(), scale, mScale.z() );
206 : 0 : }
207 : 0 : else if ( name == QLatin1String( "Z" ) )
208 : : {
209 : 0 : mOffset.set( mOffset.x(), mOffset.y(), offset );
210 : 0 : mScale.set( mScale.x(), mScale.y(), scale );
211 : 0 : }
212 : :
213 : : // store any metadata stats which are present for the attribute
214 : 0 : AttributeStatistics stats;
215 : 0 : if ( schemaObj.contains( QLatin1String( "count" ) ) )
216 : 0 : stats.count = schemaObj.value( QLatin1String( "count" ) ).toInt();
217 : 0 : if ( schemaObj.contains( QLatin1String( "minimum" ) ) )
218 : 0 : stats.minimum = schemaObj.value( QLatin1String( "minimum" ) ).toDouble();
219 : 0 : if ( schemaObj.contains( QLatin1String( "maximum" ) ) )
220 : 0 : stats.maximum = schemaObj.value( QLatin1String( "maximum" ) ).toDouble();
221 : 0 : if ( schemaObj.contains( QLatin1String( "count" ) ) )
222 : 0 : stats.mean = schemaObj.value( QLatin1String( "mean" ) ).toDouble();
223 : 0 : if ( schemaObj.contains( QLatin1String( "stddev" ) ) )
224 : 0 : stats.stDev = schemaObj.value( QLatin1String( "stddev" ) ).toDouble();
225 : 0 : if ( schemaObj.contains( QLatin1String( "variance" ) ) )
226 : 0 : stats.variance = schemaObj.value( QLatin1String( "variance" ) ).toDouble();
227 : 0 : mMetadataStats.insert( name, stats );
228 : :
229 : 0 : if ( schemaObj.contains( QLatin1String( "counts" ) ) )
230 : : {
231 : 0 : QMap< int, int > classCounts;
232 : 0 : const QJsonArray counts = schemaObj.value( QLatin1String( "counts" ) ).toArray();
233 : 0 : for ( const QJsonValue &count : counts )
234 : : {
235 : 0 : const QJsonObject countObj = count.toObject();
236 : 0 : classCounts.insert( countObj.value( QLatin1String( "value" ) ).toInt(), countObj.value( QLatin1String( "count" ) ).toInt() );
237 : 0 : }
238 : 0 : mAttributeClasses.insert( name, classCounts );
239 : 0 : }
240 : 0 : }
241 : 0 : setAttributes( attributes );
242 : :
243 : : // save mRootBounds
244 : :
245 : : // bounds (cube - octree volume)
246 : 0 : double xmin = bounds[0].toDouble();
247 : 0 : double ymin = bounds[1].toDouble();
248 : 0 : double zmin = bounds[2].toDouble();
249 : 0 : double xmax = bounds[3].toDouble();
250 : 0 : double ymax = bounds[4].toDouble();
251 : 0 : double zmax = bounds[5].toDouble();
252 : :
253 : 0 : mRootBounds = QgsPointCloudDataBounds(
254 : 0 : ( xmin - mOffset.x() ) / mScale.x(),
255 : 0 : ( ymin - mOffset.y() ) / mScale.y(),
256 : 0 : ( zmin - mOffset.z() ) / mScale.z(),
257 : 0 : ( xmax - mOffset.x() ) / mScale.x(),
258 : 0 : ( ymax - mOffset.y() ) / mScale.y(),
259 : 0 : ( zmax - mOffset.z() ) / mScale.z()
260 : : );
261 : :
262 : :
263 : : #ifdef QGIS_DEBUG
264 : : double dx = xmax - xmin, dy = ymax - ymin, dz = zmax - zmin;
265 : : QgsDebugMsgLevel( QStringLiteral( "lvl0 node size in CRS units: %1 %2 %3" ).arg( dx ).arg( dy ).arg( dz ), 2 ); // all dims should be the same
266 : : QgsDebugMsgLevel( QStringLiteral( "res at lvl0 %1" ).arg( dx / mSpan ), 2 );
267 : : QgsDebugMsgLevel( QStringLiteral( "res at lvl1 %1" ).arg( dx / mSpan / 2 ), 2 );
268 : : QgsDebugMsgLevel( QStringLiteral( "res at lvl2 %1 with node size %2" ).arg( dx / mSpan / 4 ).arg( dx / 4 ), 2 );
269 : : #endif
270 : :
271 : 0 : return true;
272 : 0 : }
273 : :
274 : 0 : QgsPointCloudBlock *QgsEptPointCloudIndex::nodeData( const IndexedPointCloudNode &n, const QgsPointCloudRequest &request )
275 : : {
276 : 0 : if ( !mHierarchy.contains( n ) )
277 : 0 : return nullptr;
278 : :
279 : 0 : if ( mDataType == QLatin1String( "binary" ) )
280 : : {
281 : 0 : QString filename = QStringLiteral( "%1/ept-data/%2.bin" ).arg( mDirectory, n.toString() );
282 : 0 : return QgsEptDecoder::decompressBinary( filename, attributes(), request.attributes() );
283 : 0 : }
284 : 0 : else if ( mDataType == QLatin1String( "zstandard" ) )
285 : : {
286 : 0 : QString filename = QStringLiteral( "%1/ept-data/%2.zst" ).arg( mDirectory, n.toString() );
287 : 0 : return QgsEptDecoder::decompressZStandard( filename, attributes(), request.attributes() );
288 : 0 : }
289 : 0 : else if ( mDataType == QLatin1String( "laszip" ) )
290 : : {
291 : 0 : QString filename = QStringLiteral( "%1/ept-data/%2.laz" ).arg( mDirectory, n.toString() );
292 : 0 : return QgsEptDecoder::decompressLaz( filename, attributes(), request.attributes() );
293 : 0 : }
294 : : else
295 : : {
296 : 0 : return nullptr; // unsupported
297 : : }
298 : 0 : }
299 : :
300 : 0 : QgsPointCloudBlockRequest *QgsEptPointCloudIndex::asyncNodeData( const IndexedPointCloudNode &n, const QgsPointCloudRequest &request )
301 : : {
302 : 0 : Q_UNUSED( n );
303 : 0 : Q_UNUSED( request );
304 : : Q_ASSERT( false );
305 : 0 : return nullptr; // unsupported
306 : : }
307 : :
308 : 0 : QgsCoordinateReferenceSystem QgsEptPointCloudIndex::crs() const
309 : : {
310 : 0 : return QgsCoordinateReferenceSystem::fromWkt( mWkt );
311 : : }
312 : :
313 : 0 : qint64 QgsEptPointCloudIndex::pointCount() const
314 : : {
315 : 0 : return mPointCount;
316 : : }
317 : :
318 : 0 : QVariant QgsEptPointCloudIndex::metadataStatistic( const QString &attribute, QgsStatisticalSummary::Statistic statistic ) const
319 : : {
320 : 0 : if ( !mMetadataStats.contains( attribute ) )
321 : 0 : return QVariant();
322 : :
323 : 0 : const AttributeStatistics &stats = mMetadataStats[ attribute ];
324 : 0 : switch ( statistic )
325 : : {
326 : : case QgsStatisticalSummary::Count:
327 : 0 : return stats.count >= 0 ? QVariant( stats.count ) : QVariant();
328 : :
329 : : case QgsStatisticalSummary::Mean:
330 : 0 : return std::isnan( stats.mean ) ? QVariant() : QVariant( stats.mean );
331 : :
332 : : case QgsStatisticalSummary::StDev:
333 : 0 : return std::isnan( stats.stDev ) ? QVariant() : QVariant( stats.stDev );
334 : :
335 : : case QgsStatisticalSummary::Min:
336 : 0 : return stats.minimum;
337 : :
338 : : case QgsStatisticalSummary::Max:
339 : 0 : return stats.maximum;
340 : :
341 : : case QgsStatisticalSummary::Range:
342 : 0 : return stats.minimum.isValid() && stats.maximum.isValid() ? QVariant( stats.maximum.toDouble() - stats.minimum.toDouble() ) : QVariant();
343 : :
344 : : case QgsStatisticalSummary::CountMissing:
345 : : case QgsStatisticalSummary::Sum:
346 : : case QgsStatisticalSummary::Median:
347 : : case QgsStatisticalSummary::StDevSample:
348 : : case QgsStatisticalSummary::Minority:
349 : : case QgsStatisticalSummary::Majority:
350 : : case QgsStatisticalSummary::Variety:
351 : : case QgsStatisticalSummary::FirstQuartile:
352 : : case QgsStatisticalSummary::ThirdQuartile:
353 : : case QgsStatisticalSummary::InterQuartileRange:
354 : : case QgsStatisticalSummary::First:
355 : : case QgsStatisticalSummary::Last:
356 : : case QgsStatisticalSummary::All:
357 : 0 : return QVariant();
358 : : }
359 : 0 : return QVariant();
360 : 0 : }
361 : :
362 : 0 : QVariantList QgsEptPointCloudIndex::metadataClasses( const QString &attribute ) const
363 : : {
364 : 0 : QVariantList classes;
365 : 0 : const QMap< int, int > values = mAttributeClasses.value( attribute );
366 : 0 : for ( auto it = values.constBegin(); it != values.constEnd(); ++it )
367 : : {
368 : 0 : classes << it.key();
369 : 0 : }
370 : 0 : return classes;
371 : 0 : }
372 : :
373 : 0 : QVariant QgsEptPointCloudIndex::metadataClassStatistic( const QString &attribute, const QVariant &value, QgsStatisticalSummary::Statistic statistic ) const
374 : : {
375 : 0 : if ( statistic != QgsStatisticalSummary::Count )
376 : 0 : return QVariant();
377 : :
378 : 0 : const QMap< int, int > values = mAttributeClasses.value( attribute );
379 : 0 : if ( !values.contains( value.toInt() ) )
380 : 0 : return QVariant();
381 : 0 : return values.value( value.toInt() );
382 : 0 : }
383 : :
384 : 0 : bool QgsEptPointCloudIndex::loadHierarchy()
385 : : {
386 : 0 : QQueue<QString> queue;
387 : 0 : queue.enqueue( QStringLiteral( "0-0-0-0" ) );
388 : 0 : while ( !queue.isEmpty() )
389 : : {
390 : 0 : const QString filename = QStringLiteral( "%1/ept-hierarchy/%2.json" ).arg( mDirectory, queue.dequeue() );
391 : 0 : QFile fH( filename );
392 : 0 : if ( !fH.open( QIODevice::ReadOnly ) )
393 : : {
394 : 0 : QgsDebugMsgLevel( QStringLiteral( "unable to read hierarchy from file %1" ).arg( filename ), 2 );
395 : 0 : return false;
396 : : }
397 : :
398 : 0 : QByteArray dataJsonH = fH.readAll();
399 : : QJsonParseError errH;
400 : 0 : const QJsonDocument docH = QJsonDocument::fromJson( dataJsonH, &errH );
401 : 0 : if ( errH.error != QJsonParseError::NoError )
402 : : {
403 : 0 : QgsDebugMsgLevel( QStringLiteral( "QJsonParseError when reading hierarchy from file %1" ).arg( filename ), 2 );
404 : 0 : return false;
405 : : }
406 : :
407 : 0 : const QJsonObject rootHObj = docH.object();
408 : 0 : for ( auto it = rootHObj.constBegin(); it != rootHObj.constEnd(); ++it )
409 : : {
410 : 0 : QString nodeIdStr = it.key();
411 : 0 : int nodePointCount = it.value().toInt();
412 : 0 : if ( nodePointCount < 0 )
413 : : {
414 : 0 : queue.enqueue( nodeIdStr );
415 : 0 : }
416 : : else
417 : : {
418 : 0 : IndexedPointCloudNode nodeId = IndexedPointCloudNode::fromString( nodeIdStr );
419 : 0 : mHierarchy[nodeId] = nodePointCount;
420 : : }
421 : 0 : }
422 : 0 : }
423 : 0 : return true;
424 : 0 : }
425 : :
426 : 0 : bool QgsEptPointCloudIndex::isValid() const
427 : : {
428 : 0 : return mIsValid;
429 : : }
430 : :
431 : : ///@endcond
|