DM_spatial_processing.cpp

//  =================================================================
//   example showing how to spatially process an ODM.
//   based on a generic neighborhood definition statistical features
//   are extracted and attached as attributes to each point
//  =================================================================
 
#include <iostream>
#include <chrono>
#include "DM/IAddInfoLayoutFactory.hpp"
#include "DM/IQueryDescriptor.hpp"
#include "DM/Datamanager.hpp"
#include "DM/Processor/IProcessorEx.hpp"
 
using namespace std;
 
// callback object computing the user-defined attributes '_pcount' and '_zmean'
class Kernel : public DM::IKernelPointEx
{
public:
  Kernel() {
  }
 
  ///returns a copy of the object (required for multi-threaded execution)
  DM::IKernelBase* threadClone() const {
    return new Kernel(*this);
  }
 
  // callback notifies about a changed leaf
  void leafChanged(const DM::IPointIndexLeaf &leaf, DM::PointIndexLeafHandle &localTree) {
    cout << "process tile with id = " << leaf.id() << endl;
  }
 
  // callback notifies that all leaf handels should be released
  void releaseLeaf() {
  }
 
  // callback for processing a point
  bool process(const DM::IPoint &source, DM::IPoint &target, const DM::IPointSet &neighbours) {
    // add neighbour count(_pcount) and z mean of neighbours(_zmean) to the current point
    double zsum = 0;
    for (unsigned n = 0; n < neighbours.sizePoint(); n++)
    {
      zsum += neighbours[n].z();
    }
 
    target.info().setInt32(0, neighbours.sizePoint());
    if (neighbours.sizePoint() > 0)
    {
      float zmean = zsum / float(neighbours.sizePoint());
      target.info().setFloat(1, zmean);
    }
 
    // return true if the point was changed otherwise false
    // this internally marks the current point leaf as changed(= needs to be written to disk)
    return true;
  }
};
 
 
void DM_spatial_processing(const char *odm)
{
  //open the odm
  DM::DatamanagerHandle dm = DM::IDatamanager::load(odm, false/*readOnly*/, true/*threadSafety*/);
 
  //if the dm wasn't opened successful exit function
  if (!dm)
  {
    cout << "Unable to open ODM '" << odm << "'" << endl;
    return;
  }
 
  // creates the appropriate attribute layout for processing
  DM::AddInfoLayoutFactoryHandle lf = DM::IAddInfoLayoutFactory::New();
  lf->addColumn(DM::ColumnType::int32, "_pcount");  // column 0
  lf->addColumn(DM::ColumnType::float_, "_zmean");  // column 1
  DM::AddInfoLayoutHandle layout = lf->getLayout();
  DM::AddInfoLayoutHandle kernelLayout = lf->getLayout();
 
  // create spatial selection
  DM::QueryDescriptorHandle query = DM::IQueryDescriptor::New("sphere(r=1)", 
      DM::IQueryDescriptor::Vocabulary::full);  // use a sphere with radius 1m as neighbourhood definition(textual description)
 
  //DM::QuerySphereHandle sphere = DM::IQuerySphere::New(/*r*/1);                 // sphere with radius 1m(class based description)
  //DM::QueryKnnHandle knn = DM::IQueryKnn::New(/*k*/12, /*dim*/3);               // 3d knn with 12 neighbours(class based description)
  //DM::QueryDescriptorHandle query = DM::IQueryDescriptor::New(sphere);          // create actual neighbourhood definition : fixed distance
  //DM::QueryDescriptorHandle query = DM::IQueryDescriptor::New(knn && sphere);   // create actual neighbourhood definition : knn with a maximum search radius
 
  cout << "Ready to append attribute '_pcount' and '_zmean' to '" << odm << "'" << endl;
 
  DM::PointIndexHandle pi = dm->getPointIndex();
  cout << "Number of tiles to process = " << pi->sizeLeaf();
 
  // create kernel object
  Kernel k;
 
  auto start = chrono::high_resolution_clock::now();
  // create processor
  // pyDM.ProcessorEx parameter : dm(datamnager), query(QueryDescriptor),
  //                              processFilter, processLayout, processLayoutReadOnly,
  //                              neighbourFilter, neighbourLayout, neighbourLayoutReadOnly
  DM::ProcessorExHandle processor = DM::IProcessorEx::New(dm, query, 
                                                          DM::FilterHandle(), layout, false, 
                                                          DM::FilterHandle(), DM::AddInfoLayoutHandle(), false);
  processor->runThreaded(k);  // perform computation
  
  auto stop = chrono::high_resolution_clock::now();
  cout << "Processing took " << chrono::duration_cast<chrono::milliseconds>(stop-start).count()/double(1000) << " [s]" << endl;
 
  // save manager object
  cout << "Save manager" << endl;
  dm->save();
 
  cout << "done" << endl;
}