From a378bcb957094f51257f5f070877ea80c3c7862d Mon Sep 17 00:00:00 2001 From: Daniel Girardeau-Montaut Date: Sat, 1 Dec 2018 22:40:04 +0100 Subject: [PATCH] Neighborhood features added (WIP) --- ContextBasedFeature.cpp | 17 ++ ContextBasedFeature.h | 6 +- DualCloudFeature.cpp | 23 +++ DualCloudFeature.h | 6 +- FeaturesInterface.cpp | 101 ++++++++++++ FeaturesInterface.h | 62 +++----- NeighborhoodFeature.cpp | 336 +++++++++++++++++++++++++++++++++++++++- NeighborhoodFeature.h | 25 +-- PointFeature.cpp | 149 ++---------------- PointFeature.h | 10 +- q3DMASCTools.cpp | 135 ++++++++++++---- 11 files changed, 638 insertions(+), 232 deletions(-) create mode 100644 FeaturesInterface.cpp diff --git a/ContextBasedFeature.cpp b/ContextBasedFeature.cpp index ea4bd08..e3d74e2 100644 --- a/ContextBasedFeature.cpp +++ b/ContextBasedFeature.cpp @@ -26,3 +26,20 @@ bool ContextBasedFeature::prepare( const CorePoints& corePoints, //TODO return false; } + +bool ContextBasedFeature::checkValidity(QString &error) const +{ + if (!Feature::checkValidity(error)) + { + return false; + } + + unsigned char cloudCount = (cloud1 ? (cloud2 ? 2 : 1) : 0); + if (cloudCount < 2) + { + error = "at least two clouds are required to compute context-based features"; + return false; + } + + return true; +} \ No newline at end of file diff --git a/ContextBasedFeature.h b/ContextBasedFeature.h index 029cdf8..27b1f16 100644 --- a/ContextBasedFeature.h +++ b/ContextBasedFeature.h @@ -77,13 +77,11 @@ namespace masc scale = p_scale; } - //! Returns the feature type + //inherited from Feature virtual Type getType() const override { return Type::ContextBasedFeature; } - //! Clones this feature virtual Feature::Shared clone() const override { return Feature::Shared(new ContextBasedFeature(*this)); } - //! Prepares the feature (compute the scalar field, etc.) virtual bool prepare(const CorePoints& corePoints, QString& error, CCLib::GenericProgressCallback* progressCb = nullptr) override; - //! Returns the descriptor for this particular feature + virtual bool checkValidity(QString &error) const override; virtual QString toString() const override { //use the default keyword + number of neighbors + the scale + the context class diff --git a/DualCloudFeature.cpp b/DualCloudFeature.cpp index d1be5b5..93ab81d 100644 --- a/DualCloudFeature.cpp +++ b/DualCloudFeature.cpp @@ -26,3 +26,26 @@ bool DualCloudFeature::prepare( const CorePoints& corePoints, //TODO return false; } + +bool DualCloudFeature::checkValidity(QString &error) const +{ + if (!Feature::checkValidity(error)) + { + return false; + } + + unsigned char cloudCount = (cloud1 ? (cloud2 ? 2 : 1) : 0); + if (cloudCount < 2) + { + error = "at least two clouds are required to compute context-based features"; + return false; + } + + if (op != NO_OPERATION) + { + error = "math operations can't be defined on dual-cloud features"; + return false; + } + + return true; +} \ No newline at end of file diff --git a/DualCloudFeature.h b/DualCloudFeature.h index 65655a0..b56fa53 100644 --- a/DualCloudFeature.h +++ b/DualCloudFeature.h @@ -65,13 +65,11 @@ namespace masc : type(p_type) {} - //! Returns the feature type + //inherited from Feature virtual Type getType() const override { return Type::DualCloudFeature; } - //! Clones this feature virtual Feature::Shared clone() const override { return Feature::Shared(new DualCloudFeature(*this)); } - //! Prepares the feature (compute the scalar field, etc.) virtual bool prepare(const CorePoints& corePoints, QString& error, CCLib::GenericProgressCallback* progressCb = nullptr) override; - //! Returns the descriptor for this particular feature + virtual bool checkValidity(QString &error) const override; virtual QString toString() const override { //use the default keyword + "_SC" + the scale diff --git a/FeaturesInterface.cpp b/FeaturesInterface.cpp new file mode 100644 index 0000000..3ded6c3 --- /dev/null +++ b/FeaturesInterface.cpp @@ -0,0 +1,101 @@ +//########################################################################## +//# # +//# CLOUDCOMPARE PLUGIN: q3DMASC # +//# # +//# This program is free software; you can redistribute it and/or modify # +//# it under the terms of the GNU General Public License as published by # +//# the Free Software Foundation; version 2 or later of the License. # +//# # +//# This program is distributed in the hope that it will be useful, # +//# but WITHOUT ANY WARRANTY; without even the implied warranty of # +//# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the # +//# GNU General Public License for more details. # +//# # +//# COPYRIGHT: Dimitri Lague / CNRS / UEB # +//# # +//########################################################################## + +#include "FeaturesInterface.h" + +//qCC_db +#include + +//system +#include + +using namespace masc; + +CCLib::ScalarField* Feature::PrepareSF(ccPointCloud* cloud, const char* resultSFName) +{ + if (!cloud || !resultSFName) + { + //invalid input parameters + assert(false); + return nullptr; + } + + CCLib::ScalarField* resultSF = nullptr; + int sfIdx = cloud->getScalarFieldIndexByName(resultSFName); + if (sfIdx >= 0) + { + resultSF = cloud->getScalarField(sfIdx); + } + else + { + ccScalarField* newSF = new ccScalarField(resultSFName); + if (!newSF->resizeSafe(cloud->size())) + { + ccLog::Warning("Not enough memory"); + newSF->release(); + return nullptr; + } + cloud->addScalarField(newSF); + + resultSF = newSF; + + } + + assert(resultSF); + resultSF->fill(NAN_VALUE); + + return resultSF; +} + +bool Feature::PerformMathOp(CCLib::ScalarField* sf1, const CCLib::ScalarField* sf2, Feature::Operation op) +{ + if (!sf1 || !sf2 || sf1->size() != sf2->size() || op == Feature::NO_OPERATION) + { + //invalid input parameters + return false; + } + + for (unsigned i = 0; i < sf1->size(); ++i) + { + ScalarType s1 = sf1->getValue(i); + ScalarType s2 = sf2->getValue(i); + ScalarType s = NAN_VALUE; + switch (op) + { + case Feature::MINUS: + s = s1 - s2; + break; + case Feature::PLUS: + s = s1 + s2; + break; + case Feature::DIVIDE: + if (std::abs(s2) > std::numeric_limits::epsilon()) + s = s1 / s2; + break; + case Feature::MULTIPLY: + s = s1 * s2; + break; + default: + assert(false); + break; + } + sf1->setValue(i, s); + } + sf1->computeMinAndMax(); + + return true; +} diff --git a/FeaturesInterface.h b/FeaturesInterface.h index cfb463b..4cad8d3 100644 --- a/FeaturesInterface.h +++ b/FeaturesInterface.h @@ -29,6 +29,11 @@ class ccPointCloud; +namespace CCLib +{ + class ScalarField; +}; + namespace masc { //! Generic feature descriptor @@ -138,6 +143,9 @@ namespace masc //! Prepares the feature (compute the scalar field, etc.) virtual bool prepare(const CorePoints& corePoints, QString& error, CCLib::GenericProgressCallback* progressCb = nullptr) = 0; + //! Finishes the feature preparation (update the scalar field, etc.) + virtual bool finish(const CorePoints& corePoints, QString& error) { /* does nothing by default*/return true; } + //! Returns whether the feature has an associated scale inline bool scaled() const { return std::isfinite(scale); } @@ -145,63 +153,35 @@ namespace masc virtual bool checkValidity(QString &error) const { unsigned char cloudCount = (cloud1 ? (cloud2 ? 2 : 1) : 0); - if (cloudCount == 0) { error = "feature has no associated cloud"; return false; } - if (scaled() && stat == NO_STAT) + if (stat != NO_STAT && getType() != Type::PointFeature) { - error = "scaled features need a STAT measure to be defined"; + error = "STAT measures can only be defined on Point features"; return false; } - if (stat != NO_STAT) + if (op != NO_OPERATION && cloudCount < 2) { - if (getType() != Type::PointFeature) - { - error = "STAT measures can only be defined on Point features"; - return false; - } - if (!scaled()) - { - error = "STAT measures need at least one scale to be defined"; - return false; - } - } - - if (op != NO_OPERATION) - { - if (!scaled()) - { - error = "math operations can't be defined on scale-less features (SC0)"; - return false; - } - if (getType() == Type::DualCloudFeature) - { - error = "math operations can't be defined on dual-cloud features"; - return false; - } - if (cloudCount < 2) - { - error = "at least two clouds are required to apply math operations"; - return false; - } - } - if (getType() == Feature::Type::DualCloudFeature || getType() == Feature::Type::ContextBasedFeature) - { - if (cloudCount < 2) - { - error = "at least two clouds are required to compute dual-cloud or context-based features"; - return false; - } + error = "at least two clouds are required to apply math operations"; + return false; } return true; } + public: //helpers + + //! Creates (or resets) a scalar field with the given name on the input core points cloud + static CCLib::ScalarField* PrepareSF(ccPointCloud* cloud, const char* resultSFName); + + //! Performs a mathematical operation between two scalar fields (they must have the same size!) + static bool PerformMathOp(CCLib::ScalarField* sf1, const CCLib::ScalarField* sf2, Operation op); + public: //members //! Scale (diameter) diff --git a/NeighborhoodFeature.cpp b/NeighborhoodFeature.cpp index b90deee..b338dcf 100644 --- a/NeighborhoodFeature.cpp +++ b/NeighborhoodFeature.cpp @@ -17,12 +17,344 @@ #include "NeighborhoodFeature.h" +//CCLib +#include +#include + using namespace masc; +bool NeighborhoodFeature::checkValidity(QString &error) const +{ + if (!Feature::checkValidity(error)) + { + return false; + } + + if (stat != Feature::NO_STAT) + { + error = "Neighborhood features shouldn't be associated to a STAT measure"; + return false; + } + + if (cloud2 && op == NO_OPERATION) + { + error = "Feature has a second cloud associated but no MATH operation is defined"; + return false; + } + + return true; +} + bool NeighborhoodFeature::prepare( const CorePoints& corePoints, QString& error, CCLib::GenericProgressCallback* progressCb/*=nullptr*/) { - //TODO - return false; + if (!cloud1 || !corePoints.cloud) + { + //invalid input + assert(false); + error = "internal error (no input core points)"; + return false; + } + + if (!checkValidity(error)) + { + assert(false); + return false; + } + + //build the final SF name + QString resultSFName = ToString(type) + "_" + cloud1Label; + if (cloud2) + { + //include the math operation as well if necessary! + resultSFName += "_" + Feature::OpToString(op) + "_" + cloud2Label; + } + resultSFName += "@" + QString::number(scale); + + //and the scalar field + assert(!sf1); + sf1 = PrepareSF(corePoints.cloud, qPrintable(resultSFName)); + if (!sf1) + { + error = QString("Failed to prepare scalar %1 @ scale %2").arg(cloud1Label).arg(scale); + return false; + } + sourceName = sf1->getName(); + + if (cloud2 && op != Feature::NO_OPERATION) + { + QString resultSFName2 = ToString(type) + "_" + cloud2Label + "@" + QString::number(scale); + keepSF2 = (corePoints.cloud->getScalarFieldIndexByName(qPrintable(resultSFName2)) >= 0); //we remember that the scalar field was already existing! + + assert(!sf2); + sf2 = PrepareSF(corePoints.cloud, qPrintable(resultSFName2)); + if (!sf2) + { + error = QString("Failed to prepare scalar field for %1 @ scale %2").arg(cloud2Label).arg(scale); + return false; + } + } + + return true; +} + +bool NeighborhoodFeature::finish(const CorePoints& corePoints, QString& error) +{ + if (!corePoints.cloud) + { + //invalid input + assert(false); + error = "internal error (no input core points)"; + return false; + } + + bool success = true; + + if (sf1) + { + sf1->computeMinAndMax(); + + //update display + //if (corePoints.cloud->getDisplay()) + { + int sfIndex1 = corePoints.cloud->getScalarFieldIndexByName(sf1->getName()); + corePoints.cloud->setCurrentDisplayedScalarField(sfIndex1); + //corePoints.cloud->getDisplay()->redraw(); + //QCoreApplication::processEvents(); + } + } + + if (sf2) + { + //now perform the math operation + if (op != Feature::NO_OPERATION) + { + if (!PerformMathOp(sf1, sf2, op)) + { + error = "Failed to perform the MATH operation"; + success = false; + } + } + + if (keepSF2) + { + sf2->computeMinAndMax(); + } + else + { + int sfIndex2 = corePoints.cloud->getScalarFieldIndexByName(sf2->getName()); + if (sfIndex2 >= 0) + { + corePoints.cloud->deleteScalarField(sfIndex2); + } + else + { + assert(false); + sf2->release(); + } + sf2 = nullptr; + } + } + + return success; +} + +QString NeighborhoodFeature::toString() const +{ + //use the default keyword + the scale + QString description = ToString(type) + "_SC" + QString::number(scale); + + description += "_" + cloud1Label; + + if (cloud2 && !cloud2Label.isEmpty()) + { + description += "_" + cloud2Label; + + if (op != NO_OPERATION) + { + description += "_" + OpToString(op); + } + } + + return description; +} + +bool NeighborhoodFeature::computeValue(CCLib::DgmOctree::NeighboursSet& pointsInNeighbourhood, const CCVector3& queryPoint, double& outputValue) const +{ + outputValue = std::numeric_limits::quiet_NaN(); + + size_t kNN = pointsInNeighbourhood.size(); + if (kNN == 0) + { + assert(false); + return false; + } + + switch (type) + { + //features relying on the PCA + case PCA1: + case PCA2: + case SPHER: + case LINEA: + case PLANA: + case FOM: + case LINEF: + case ORIENF: + if (kNN >= 6) + { + CCLib::DgmOctreeReferenceCloud neighboursCloud(&pointsInNeighbourhood, static_cast(kNN)); + CCLib::Neighbourhood Z(&neighboursCloud); + CCLib::SquareMatrixd eigVectors; + std::vector eigValues; + if (Jacobi::ComputeEigenValuesAndVectors(Z.computeCovarianceMatrix(), eigVectors, eigValues, true)) + { + Jacobi::SortEigenValuesAndVectors(eigVectors, eigValues); //decreasing order of their associated eigenvalues + switch (type) + { + case PCA1: + outputValue = eigValues[0] / (eigValues[0] + eigValues[1] + eigValues[2]); + break; + case PCA2: + outputValue = eigValues[1] / (eigValues[0] + eigValues[1] + eigValues[2]); + break; + case SPHER: + if (std::abs(eigValues[0]) > std::numeric_limits::epsilon()) + outputValue = eigValues[2] / eigValues[0]; + break; + case LINEA: + if (std::abs(eigValues[0]) > std::numeric_limits::epsilon()) + outputValue = (eigValues[0] - eigValues[1]) / eigValues[0]; + break; + case PLANA: + if (std::abs(eigValues[0]) > std::numeric_limits::epsilon()) + outputValue = (eigValues[1] - eigValues[2]) / eigValues[0]; + break; + case FOM: + { + double m1 = 0.0, m2 = 0.0; + CCVector3 e2(eigVectors.m_values[0][1], eigVectors.m_values[1][1], eigVectors.m_values[2][1]); + for (size_t i = 0; i < kNN; ++i) + { + double dotProd = (*pointsInNeighbourhood[i].point - queryPoint).dot(e2); + m1 += dotProd; + m2 += dotProd * dotProd; + } + outputValue = (m1 * m1) / m2; + } + case LINEF: + case ORIENF: + //can't compute these values yet! + break; + default: + //impossible + assert(false); + break; + } + } + else + { + return false; + } + } + break; + + case DipAng: + case DipDir: + if (kNN >= 3) + { + CCLib::DgmOctreeReferenceCloud neighboursCloud(&pointsInNeighbourhood, static_cast(kNN)); + CCLib::Neighbourhood Z(&neighboursCloud); + const CCVector3* N = Z.getLSPlaneNormal(); + if (N) + { + //force +Z + CCVector3 Np = (N->z < 0 ? -PC_ONE * *N : *N); + PointCoordinateType dip_deg, dipDir_deg; + ccNormalVectors::ConvertNormalToDipAndDipDir(Np, dip_deg, dipDir_deg); + outputValue = (type == DipAng ? dip_deg : dipDir_deg); + } + else + { + return false; + } + } + break; + + case ROUGH: + + case NBPTS: + outputValue = static_cast(kNN); + break; + + case CURV: + + case ZRANGE: + case Zmax: + case Zmin: + if (kNN >= 2) + { + PointCoordinateType minZ, maxZ; + minZ = maxZ = pointsInNeighbourhood[0].point->z; + for (size_t i = 1; i < kNN; ++i) + { + if (minZ < pointsInNeighbourhood[i].point->z) + minZ = pointsInNeighbourhood[i].point->z; + else if (maxZ > pointsInNeighbourhood[i].point->z) + maxZ = pointsInNeighbourhood[i].point->z; + } + + if (type == ZRANGE) + { + outputValue = maxZ - minZ; + } + else if (type == Zmax) + { + outputValue = maxZ - queryPoint.z; + } + else if (type == Zmax) + { + outputValue = queryPoint.z - minZ; + } + else + { + //impossible + assert(false); + } + } + + case ANISO: + if (kNN >= 3) + { + CCLib::DgmOctreeReferenceCloud neighboursCloud(&pointsInNeighbourhood, static_cast(kNN)); + CCLib::Neighbourhood Z(&neighboursCloud); + const CCVector3* G = Z.getGravityCenter(); + if (G) + { + double r = sqrt(pointsInNeighbourhood.back().squareDistd); + if (r > std::numeric_limits::epsilon()) + { + double d = (queryPoint - *G).normd(); + //Ratio of distance to center of mass and radius of sphere + outputValue = d / r; + } + } + else + { + return false; + } + } + break; + + default: + { + ccLog::Warning("Unhandled STAT measure"); + assert(false); + return false; + } + + } + + return true; } diff --git a/NeighborhoodFeature.h b/NeighborhoodFeature.h index df9f536..39599ba 100644 --- a/NeighborhoodFeature.h +++ b/NeighborhoodFeature.h @@ -27,6 +27,8 @@ namespace masc { public: //NeighborhoodFeatureType + typedef QSharedPointer Shared; + enum NeighborhoodFeatureType { Invalid = 0 @@ -143,21 +145,22 @@ namespace masc //! Default constructor NeighborhoodFeature(NeighborhoodFeatureType p_type) : type(p_type) + , sf1(nullptr) + , sf2(nullptr) + , keepSF2(false) { } - //! Returns the feature type + //inherited from Feature virtual Type getType() const override { return Type::NeighborhoodFeature; } - //! Clones this feature virtual Feature::Shared clone() const override { return Feature::Shared(new NeighborhoodFeature(*this)); } - //! Prepares the feature (compute the scalar field, etc.) virtual bool prepare(const CorePoints& corePoints, QString& error, CCLib::GenericProgressCallback* progressCb = nullptr) override; - //! Returns the descriptor for this particular feature - virtual QString toString() const override - { - //use the default keyword + the scale - return ToString(type) + "_SC" + QString::number(scale); - } + virtual bool finish(const CorePoints& corePoints, QString& error) override; + virtual bool checkValidity(QString &error) const override; + virtual QString toString() const override; + + //! Compute the feature values on a set of points + bool computeValue(CCLib::DgmOctree::NeighboursSet& pointsInNeighbourhood, const CCVector3& queryPoint, double& outputValue) const; public: //members @@ -165,5 +168,9 @@ namespace masc /** \warning different from the feature type **/ NeighborhoodFeatureType type; + + //! Feature values + CCLib::ScalarField *sf1, *sf2; + bool keepSF2; }; } \ No newline at end of file diff --git a/PointFeature.cpp b/PointFeature.cpp index dc6eea0..8b1099c 100644 --- a/PointFeature.cpp +++ b/PointFeature.cpp @@ -54,6 +54,18 @@ bool PointFeature::checkValidity(QString &error) const assert(cloud1); + if (scaled() && stat == NO_STAT) + { + error = "scaled point features need a STAT measure to be defined"; + return false; + } + + if (op != NO_OPERATION && !scaled()) + { + error = "math operations can't be defined on scale-less point features (SC0)"; + return false; + } + switch (type) { case PointFeature::Intensity: @@ -439,66 +451,6 @@ static bool ExtractStatFromSF( const CCVector3& queryPoint, return true; } -static bool PrepareOctree(ccPointCloud* sourceCloud, CCLib::GenericProgressCallback* progressCb = nullptr) -{ - if (!sourceCloud) - { - //invalid input parameters - assert(false); - return false; - } - - ccOctree::Shared octree = sourceCloud->getOctree(); - if (!octree) - { - ccLog::Print(QString("Computing octree of cloud %1 (%2 points)").arg(sourceCloud->getName()).arg(sourceCloud->size())); - octree = sourceCloud->computeOctree(progressCb); - if (!octree) - { - ccLog::Warning("Failed to compute octree"); - return nullptr; - } - } - - return true; -} - -static CCLib::ScalarField* PrepareSF(const CorePoints& corePoints, const char* resultSFName) -{ - if (!corePoints.cloud || !resultSFName) - { - //invalid input parameters - assert(false); - return nullptr; - } - - CCLib::ScalarField* resultSF = nullptr; - int sfIdx = corePoints.cloud->getScalarFieldIndexByName(resultSFName); - if (sfIdx >= 0) - { - resultSF = corePoints.cloud->getScalarField(sfIdx); - } - else - { - ccScalarField* newSF = new ccScalarField(resultSFName); - if (!newSF->resizeSafe(corePoints.cloud->size())) - { - ccLog::Warning("Not enough memory"); - newSF->release(); - return nullptr; - } - corePoints.cloud->addScalarField(newSF); - - resultSF = newSF; - - } - - assert(resultSF); - resultSF->fill(NAN_VALUE); - - return resultSF; -} - static CCLib::ScalarField* ExtractStat( const CorePoints& corePoints, ccPointCloud* sourceCloud, const IScalarFieldWrapper* sourceField, @@ -596,45 +548,6 @@ static CCLib::ScalarField* ExtractStat( const CorePoints& corePoints, return resultSF; } -static bool PerformMathOp(CCLib::ScalarField* sf1, const CCLib::ScalarField* sf2, Feature::Operation op) -{ - if (!sf1 || !sf2 || sf1->size() != sf2->size() || op == Feature::NO_OPERATION) - { - //invalid input parameters - return false; - } - - for (unsigned i = 0; i < sf1->size(); ++i) - { - ScalarType s1 = sf1->getValue(i); - ScalarType s2 = sf2->getValue(i); - ScalarType s = NAN_VALUE; - switch (op) - { - case Feature::MINUS: - s = s1 - s2; - break; - case Feature::PLUS: - s = s1 + s2; - break; - case Feature::DIVIDE: - if (std::abs(s2) > std::numeric_limits::epsilon()) - s = s1 / s2; - break; - case Feature::MULTIPLY: - s = s1 * s2; - break; - default: - assert(false); - break; - } - sf1->setValue(i, s); - } - sf1->computeMinAndMax(); - - return true; -} - bool PointFeature::prepare( const CorePoints& corePoints, QString& error, CCLib::GenericProgressCallback* progressCb/*=nullptr*/) @@ -683,6 +596,7 @@ bool PointFeature::prepare( const CorePoints& corePoints, { assert(false); ccLog::Warning("Feature has a second cloud associated but no MATH operation is defined"); + return false; } } @@ -695,20 +609,11 @@ bool PointFeature::prepare( const CorePoints& corePoints, } resultSFName += "@" + QString::number(scale); - //prepare the octree - //if (!PrepareOctree(cloud1, progressCb)) - //{ - // error = "Failed to compute octree (not enough memory?)"; - // return false; - //} - //and the scalar field assert(!statSF1); - statSF1 = PrepareSF(corePoints, qPrintable(resultSFName)); - //CCLib::ScalarField* statSF1 = ExtractStat(corePoints, cloud1, field1.data(), scale, stat, qPrintable(resultSFName), progressCb); + statSF1 = PrepareSF(corePoints.cloud, qPrintable(resultSFName)); if (!statSF1) { - //error = QString("Failed to extract stat. from field '%1' @ scale %2").arg(field1->getName()).arg(scale); error = QString("Failed to prepare scalar field for field '%1' @ scale %2").arg(field1->getName()).arg(scale); return false; } @@ -719,35 +624,13 @@ bool PointFeature::prepare( const CorePoints& corePoints, QString resultSFName2 = cloud2Label + "." + field2->getName() + QString("_") + Feature::StatToString(stat) + "@" + QString::number(scale); keepStatSF2 = (corePoints.cloud->getScalarFieldIndexByName(qPrintable(resultSFName2)) >= 0); //we remember that the scalar field was already existing! - //prepare the octree - //if (!PrepareOctree(cloud2, progressCb)) - //{ - // error = "Failed to compute octree (not enough memory?)"; - // return false; - //} - assert(!statSF2); - statSF2 = PrepareSF(corePoints, qPrintable(resultSFName2)); - //statSF2 = ExtractStat(corePoints, cloud2, field2.data(), scale, stat, qPrintable(resultSFName2), progressCb); + statSF2 = PrepareSF(corePoints.cloud, qPrintable(resultSFName2)); if (!statSF2) { - error = QString("Failed to extract stat. from field '%1' @ scale %2").arg(field2->getName()).arg(scale); + error = QString("Failed to prepare scalar field for field '%1' @ scale %2").arg(field2->getName()).arg(scale); return false; } - - //now perform the math operation - //if (!PerformMathOp(statSF1, statSF2, op)) - //{ - // error = "Failed to perform the MATH operation"; - // return false; - //} - - //if (sfIndex2 < 0) - //{ - // //release some memory - // sfIndex2 = corePoints.cloud->getScalarFieldIndexByName(qPrintable(resultSFName2)); - // corePoints.cloud->deleteScalarField(sfIndex2); - //} } return true; diff --git a/PointFeature.h b/PointFeature.h index 9734b9b..5c9492d 100644 --- a/PointFeature.h +++ b/PointFeature.h @@ -184,20 +184,14 @@ namespace masc } } - //! Returns the feature type + //inherited from Feature virtual Type getType() const override { return Type::PointFeature; } - //! Clones this feature virtual Feature::Shared clone() const override { return Feature::Shared(new PointFeature(*this)); } - //! Prepares the feature (compute the scalar field, etc.) virtual bool prepare(const CorePoints& corePoints, QString& error, CCLib::GenericProgressCallback* progressCb = nullptr) override; - //! Checks the feature definition validity + virtual bool finish(const CorePoints& corePoints, QString& error) override; virtual bool checkValidity(QString &error) const override; - //! Returns the descriptor for this particular feature virtual QString toString() const override; - //! Finishes the feature preparation (update the scalar field, etc.) - bool finish(const CorePoints& corePoints, QString& error); - //! Compute the associated 'stat' on a set of points (and with a given field) bool computeStat(const CCLib::DgmOctree::NeighboursSet& pointsInNeighbourhood, const QSharedPointer& sourceField, double& outputValue) const; diff --git a/q3DMASCTools.cpp b/q3DMASCTools.cpp index 52bc6d3..a78f008 100644 --- a/q3DMASCTools.cpp +++ b/q3DMASCTools.cpp @@ -751,6 +751,8 @@ bool Tools::LoadTrainingFile( QString filename, { loadedClouds.push_back(it.value()); } + + return true; } else { @@ -796,8 +798,9 @@ CCLib::ScalarField* Tools::RetrieveSF(const ccPointCloud* cloud, const QString& struct FeaturesAndScales { std::vector scales; - std::vector pointFeatures; - std::vector neighborhoodFeatures; + size_t featureCount = 0; + QMap > pointFeaturesPerScale; + QMap > neighborhoodFeaturesPerScale; }; bool Tools::PrepareFeatures(const CorePoints& corePoints, Feature::Set& features, QString& error, CCLib::GenericProgressCallback* progressCb/*=nullptr*/) @@ -811,7 +814,7 @@ bool Tools::PrepareFeatures(const CorePoints& corePoints, Feature::Set& features //gather all the scales that need to be extracted QMap cloudsWithScaledFeatures; - + //and prepare the features (scalar fields, etc.) at the same time for (const Feature::Shared& feature : features) { QString errorMessage("invalid pointer"); @@ -828,30 +831,73 @@ bool Tools::PrepareFeatures(const CorePoints& corePoints, Feature::Set& features return false; } - if (feature->getType() == Feature::Type::PointFeature && feature->scaled()) + if (feature->scaled()) { try { - //build the scaled feature list attached to the first cloud - if (feature->cloud1) + switch (feature->getType()) { - FeaturesAndScales& fas = cloudsWithScaledFeatures[feature->cloud1]; - fas.pointFeatures.push_back(qSharedPointerCast(feature)); - if (std::find(fas.scales.begin(), fas.scales.end(), feature->scale) == fas.scales.end()) + //Point features + case Feature::Type::PointFeature: + { + //build the scaled feature list attached to the first cloud + if (feature->cloud1) { - fas.scales.push_back(feature->scale); + FeaturesAndScales& fas = cloudsWithScaledFeatures[feature->cloud1]; + fas.pointFeaturesPerScale[feature->scale].push_back(qSharedPointerCast(feature)); + ++fas.featureCount; + if (std::find(fas.scales.begin(), fas.scales.end(), feature->scale) == fas.scales.end()) + { + fas.scales.push_back(feature->scale); + } + } + + //build the scaled feature list attached to the second cloud (if any) + if (feature->cloud2 && feature->cloud2 != feature->cloud1 && feature->op != Feature::NO_OPERATION) + { + FeaturesAndScales& fas = cloudsWithScaledFeatures[feature->cloud2]; + ++fas.featureCount; + fas.pointFeaturesPerScale[feature->scale].push_back(qSharedPointerCast(feature)); + if (std::find(fas.scales.begin(), fas.scales.end(), feature->scale) == fas.scales.end()) + { + fas.scales.push_back(feature->scale); + } } } + break; - //build the scaled feature list attached to the second cloud (if any) - if (feature->cloud2 && feature->cloud2 != feature->cloud1 && feature->op != Feature::NO_OPERATION) + //Point features + case Feature::Type::NeighborhoodFeature: { - FeaturesAndScales& fas = cloudsWithScaledFeatures[feature->cloud2]; - fas.pointFeatures.push_back(qSharedPointerCast(feature)); - if (std::find(fas.scales.begin(), fas.scales.end(), feature->scale) == fas.scales.end()) + //build the scaled feature list attached to the first cloud + if (feature->cloud1) { - fas.scales.push_back(feature->scale); + FeaturesAndScales& fas = cloudsWithScaledFeatures[feature->cloud1]; + fas.neighborhoodFeaturesPerScale[feature->scale].push_back(qSharedPointerCast(feature)); + ++fas.featureCount; + if (std::find(fas.scales.begin(), fas.scales.end(), feature->scale) == fas.scales.end()) + { + fas.scales.push_back(feature->scale); + } } + + //build the scaled feature list attached to the second cloud (if any) + if (feature->cloud2 && feature->cloud2 != feature->cloud1 && feature->op != Feature::NO_OPERATION) + { + FeaturesAndScales& fas = cloudsWithScaledFeatures[feature->cloud2]; + fas.neighborhoodFeaturesPerScale[feature->scale].push_back(qSharedPointerCast(feature)); + ++fas.featureCount; + if (std::find(fas.scales.begin(), fas.scales.end(), feature->scale) == fas.scales.end()) + { + fas.scales.push_back(feature->scale); + } + } + } + break; + + default: + assert(false); + break; } } catch (const std::bad_alloc&) @@ -860,6 +906,7 @@ bool Tools::PrepareFeatures(const CorePoints& corePoints, Feature::Set& features return false; } } + } bool success = true; @@ -894,8 +941,7 @@ bool Tools::PrepareFeatures(const CorePoints& corePoints, Feature::Set& features unsigned char octreeLevel = octree->findBestLevelForAGivenNeighbourhoodSizeExtraction(largestRadius); unsigned pointCount = corePoints.size(); - size_t featureCount = fas.pointFeatures.size() + fas.neighborhoodFeatures.size(); - QString logMessage = QString("Computing %1 features on cloud %2\n(core points: %3)").arg(featureCount).arg(sourceCloud->getName()).arg(pointCount); + QString logMessage = QString("Computing %1 features on cloud %2\n(core points: %3)").arg(fas.featureCount).arg(sourceCloud->getName()).arg(pointCount); if (progressCb) { progressCb->setMethodTitle("Compute features"); @@ -905,8 +951,10 @@ bool Tools::PrepareFeatures(const CorePoints& corePoints, Feature::Set& features CCLib::NormalizedProgress nProgress(progressCb, pointCount); QMutex mutex; +#ifndef _DEBUG #if defined(_OPENMP) #pragma omp parallel for +#endif #endif for (int i = 0; i < static_cast(pointCount); ++i) { @@ -953,14 +1001,9 @@ bool Tools::PrepareFeatures(const CorePoints& corePoints, Feature::Set& features nNSS.pointsInNeighbourhood.resize(kNN); } - for (PointFeature::Shared& feature : fas.pointFeatures) + //Point features + for (PointFeature::Shared& feature : fas.pointFeaturesPerScale[fas.scales[scaleIndex]]) { - if (feature->scale != fas.scales[scaleIndex]) - { - //we use the current neighborhood only for the features with the corresponding scales! - continue; - } - if (feature->cloud1 == sourceCloud && feature->statSF1 && feature->field1) { double outputValue = 0; @@ -991,6 +1034,39 @@ bool Tools::PrepareFeatures(const CorePoints& corePoints, Feature::Set& features } } + //Neighborhhod features + for (NeighborhoodFeature::Shared& feature : fas.neighborhoodFeaturesPerScale[fas.scales[scaleIndex]]) + { + if (feature->cloud1 == sourceCloud && feature->sf1) + { + double outputValue = 0; + if (!feature->computeValue(nNSS.pointsInNeighbourhood, nNSS.queryPoint, outputValue)) + { + //an error occurred + success = false; + break; + } + + ScalarType v1 = static_cast(outputValue); + feature->sf1->setValue(i, v1); + } + + if (feature->cloud2 == sourceCloud && feature->sf2) + { + assert(feature->op != Feature::NO_OPERATION); + double outputValue = 0; + if (!feature->computeValue(nNSS.pointsInNeighbourhood, nNSS.queryPoint, outputValue)) + { + //an error occurred + success = false; + break; + } + + ScalarType v2 = static_cast(outputValue); + feature->sf2->setValue(i, v2); + } + } + if (!success) { break; @@ -1019,13 +1095,10 @@ bool Tools::PrepareFeatures(const CorePoints& corePoints, Feature::Set& features for (const Feature::Shared& feature : features) { - //we have to 'finish' the process for Point features - if (feature->getType() == Feature::Type::PointFeature && feature->scaled()) + //we have to 'finish' the process for scaled features + if (feature->scaled() && !feature->finish(corePoints, error)) { - if (!qSharedPointerCast(feature)->finish(corePoints, error)) - { - return false; - } + return false; } }