qgsalgorithmpointstopaths.cpp

Go to the documentation of this file.

/***************************************************************************
                         qgsalgorithmdpointstopaths.cpp
                         ---------------------
    begin                : November 2020
    copyright            : (C) 2020 by Stefanos Natsis
    email                : uclaros at gmail dot com
 ***************************************************************************/
 
/***************************************************************************
 *                                                                         *
 *   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; either version 2 of the License, or     *
 *   (at your option) any later version.                                   *
 *                                                                         *
 ***************************************************************************/
 
#include "qgsalgorithmpointstopaths.h"
 
#include "qgsdistancearea.h"
#include "qgsexpressionnodeimpl.h"
#include "qgsmultipoint.h"
#include "qgsvectorlayer.h"
 
#include <QCollator>
#include <QString>
#include <QTextStream>
 
using namespace Qt::StringLiterals;

 
QString QgsPointsToPathsAlgorithm::name() const
{
  return u"pointstopath"_s;
}
 
QString QgsPointsToPathsAlgorithm::displayName() const
{
  return QObject::tr( "Points to path" );
}
 
QString QgsPointsToPathsAlgorithm::shortHelpString() const
{
  return QObject::tr(
    "This algorithm takes a point layer and connects its features to create a new line layer.\n\n"
    "An attribute or expression may be specified to define the order the points should be connected. "
    "If no order expression is specified, the feature ID is used.\n\n"
    "A natural sort can be used when sorting by a string attribute "
    "or expression (ie. place 'a9' before 'a10').\n\n"
    "An attribute or expression can be selected to group points having the same value into the same resulting line."
  );
}
 
QString QgsPointsToPathsAlgorithm::shortDescription() const
{
  return QObject::tr( "Takes a point layer and connects its features to create a new line layer." );
}
 
Qgis::ProcessingAlgorithmDocumentationFlags QgsPointsToPathsAlgorithm::documentationFlags() const
{
  return Qgis::ProcessingAlgorithmDocumentationFlag::RespectsEllipsoid;
}
 
QStringList QgsPointsToPathsAlgorithm::tags() const
{
  return QObject::tr( "create,lines,points,connect,convert,join,path" ).split( ',' );
}
 
QString QgsPointsToPathsAlgorithm::group() const
{
  return QObject::tr( "Vector creation" );
}
 
QString QgsPointsToPathsAlgorithm::groupId() const
{
  return u"vectorcreation"_s;
}
 
void QgsPointsToPathsAlgorithm::initAlgorithm( const QVariantMap & )
{
  addParameter( std::make_unique<QgsProcessingParameterFeatureSource>( u"INPUT"_s, QObject::tr( "Input layer" ), QList<int>() << static_cast<int>( Qgis::ProcessingSourceType::VectorPoint ) ) );
  addParameter( std::make_unique<QgsProcessingParameterBoolean>( u"CLOSE_PATH"_s, QObject::tr( "Create closed paths" ), false ) );
  addParameter( std::make_unique<QgsProcessingParameterExpression>( u"ORDER_EXPRESSION"_s, QObject::tr( "Order expression" ), QVariant(), u"INPUT"_s, true ) );
  addParameter( std::make_unique<QgsProcessingParameterBoolean>( u"NATURAL_SORT"_s, QObject::tr( "Sort text containing numbers naturally" ), false ) );
  addParameter( std::make_unique<QgsProcessingParameterExpression>( u"GROUP_EXPRESSION"_s, QObject::tr( "Path group expression" ), QVariant(), u"INPUT"_s, true ) );
  addParameter( std::make_unique<QgsProcessingParameterFeatureSink>( u"OUTPUT"_s, QObject::tr( "Paths" ), Qgis::ProcessingSourceType::VectorLine ) );
  // TODO QGIS 5: remove parameter. move logic to separate algorithm if needed.
  addParameter( std::make_unique<QgsProcessingParameterFolderDestination>( u"OUTPUT_TEXT_DIR"_s, QObject::tr( "Directory for text output" ), QVariant(), true, false ) );
  addOutput( std::make_unique<QgsProcessingOutputNumber>( u"NUM_PATHS"_s, QObject::tr( "Number of paths" ) ) );
 
  // backwards compatibility parameters
  // TODO QGIS 5: remove compatibility parameters and their logic
  auto orderField = std::make_unique<QgsProcessingParameterField>( u"ORDER_FIELD"_s, QObject::tr( "Order field" ), QVariant(), QString(), Qgis::ProcessingFieldParameterDataType::Any, false, true );
  orderField->setFlags( orderField->flags() | Qgis::ProcessingParameterFlag::Hidden );
  addParameter( std::move( orderField ) );
 
  auto groupField = std::make_unique<QgsProcessingParameterField>( u"GROUP_FIELD"_s, QObject::tr( "Group field" ), QVariant(), u"INPUT"_s, Qgis::ProcessingFieldParameterDataType::Any, false, true );
  groupField->setFlags( groupField->flags() | Qgis::ProcessingParameterFlag::Hidden );
  addParameter( std::move( groupField ) );
 
  auto dateFormat = std::make_unique<QgsProcessingParameterString>( u"DATE_FORMAT"_s, QObject::tr( "Date format (if order field is DateTime)" ), QVariant(), false, true );
  dateFormat->setFlags( dateFormat->flags() | Qgis::ProcessingParameterFlag::Hidden );
  addParameter( std::move( dateFormat ) );
}
 
QgsPointsToPathsAlgorithm *QgsPointsToPathsAlgorithm::createInstance() const
{
  return new QgsPointsToPathsAlgorithm();
}
 
QVariantMap QgsPointsToPathsAlgorithm::processAlgorithm( const QVariantMap &parameters, QgsProcessingContext &context, QgsProcessingFeedback *feedback )
{
  std::unique_ptr<QgsProcessingFeatureSource> source( parameterAsSource( parameters, u"INPUT"_s, context ) );
  if ( !source )
    throw QgsProcessingException( invalidSourceError( parameters, u"INPUT"_s ) );
 
  const bool closePaths = parameterAsBool( parameters, u"CLOSE_PATH"_s, context );
 
  QString orderExpressionString = parameterAsString( parameters, u"ORDER_EXPRESSION"_s, context );
  const QString orderFieldString = parameterAsString( parameters, u"ORDER_FIELD"_s, context );
  if ( !orderFieldString.isEmpty() )
  {
    // this is a backwards compatibility parameter
    orderExpressionString = QgsExpression::quotedColumnRef( orderFieldString );
 
    QString dateFormat = parameterAsString( parameters, u"DATE_FORMAT"_s, context );
    if ( !dateFormat.isEmpty() )
    {
      QVector<QPair<QString, QString>> codeMap;
      codeMap
        << QPair<QString, QString>( "%%", "%" )
        << QPair<QString, QString>( "%a", "ddd" )
        << QPair<QString, QString>( "%A", "dddd" )
        << QPair<QString, QString>( "%w", "" ) //day of the week 0-6
        << QPair<QString, QString>( "%d", "dd" )
        << QPair<QString, QString>( "%b", "MMM" )
        << QPair<QString, QString>( "%B", "MMMM" )
        << QPair<QString, QString>( "%m", "MM" )
        << QPair<QString, QString>( "%y", "yy" )
        << QPair<QString, QString>( "%Y", "yyyy" )
        << QPair<QString, QString>( "%H", "hh" )
        << QPair<QString, QString>( "%I", "hh" ) // 12 hour
        << QPair<QString, QString>( "%p", "AP" )
        << QPair<QString, QString>( "%M", "mm" )
        << QPair<QString, QString>( "%S", "ss" )
        << QPair<QString, QString>( "%f", "zzz" ) // milliseconds instead of microseconds
        << QPair<QString, QString>( "%z", "" )    // utc offset
        << QPair<QString, QString>( "%Z", "" )    // timezone name
        << QPair<QString, QString>( "%j", "" )    // day of the year
        << QPair<QString, QString>( "%U", "" )    // week number of the year sunday based
        << QPair<QString, QString>( "%W", "" )    // week number of the year monday based
        << QPair<QString, QString>( "%c", "" )    // full datetime
        << QPair<QString, QString>( "%x", "" )    // full date
        << QPair<QString, QString>( "%X", "" )    // full time
        << QPair<QString, QString>( "%G", "yyyy" )
        << QPair<QString, QString>( "%u", "" )  // day of the week 1-7
        << QPair<QString, QString>( "%V", "" ); // week number
      for ( const auto &pair : std::as_const( codeMap ) )
      {
        dateFormat.replace( pair.first, pair.second );
      }
      orderExpressionString = QString( "to_datetime(%1, '%2')" ).arg( orderExpressionString, dateFormat );
    }
  }
  else if ( orderExpressionString.isEmpty() )
  {
    // If no order expression is given, default to the fid
    orderExpressionString = QString( "$id" );
  }
  QgsExpressionContext expressionContext = createExpressionContext( parameters, context, source.get() );
  QgsExpression orderExpression = QgsExpression( orderExpressionString );
  if ( orderExpression.hasParserError() )
    throw QgsProcessingException( orderExpression.parserErrorString() );
 
  QStringList requiredFields = QStringList( orderExpression.referencedColumns().values() );
  orderExpression.prepare( &expressionContext );
 
  QMetaType::Type orderFieldType = QMetaType::Type::QString;
  if ( orderExpression.isField() )
  {
    const QString orderField = qgis::down_cast<const QgsExpressionNodeColumnRef *>( orderExpression.rootNode() )->name();
    const int orderFieldIndex = source->fields().lookupField( orderField );
    if ( orderFieldIndex < 0 )
    {
      throw QgsProcessingException( QObject::tr( "Order by field %1 does not exist in input layer." ).arg( orderField ) );
    }
    orderFieldType = source->fields().field( orderFieldIndex ).type();
  }
 
  QString groupExpressionString = parameterAsString( parameters, u"GROUP_EXPRESSION"_s, context );
  // handle backwards compatibility parameter GROUP_FIELD
  const QString groupFieldString = parameterAsString( parameters, u"GROUP_FIELD"_s, context );
  if ( !groupFieldString.isEmpty() )
    groupExpressionString = QgsExpression::quotedColumnRef( groupFieldString );
 
  QgsExpression groupExpression = groupExpressionString.isEmpty() ? QgsExpression( QString( "true" ) ) : QgsExpression( groupExpressionString );
  if ( groupExpression.hasParserError() )
    throw QgsProcessingException( groupExpression.parserErrorString() );
 
  QgsFields outputFields = QgsFields();
  if ( !groupExpressionString.isEmpty() )
  {
    requiredFields.append( groupExpression.referencedColumns().values() );
    if ( groupExpression.isField() )
    {
      const QString groupField = qgis::down_cast<const QgsExpressionNodeColumnRef *>( groupExpression.rootNode() )->name();
      const int groupFieldIndex = source->fields().lookupField( groupField );
      if ( groupFieldIndex < 0 )
      {
        throw QgsProcessingException( QObject::tr( "Group field %1 does not exist in input layer." ).arg( groupField ) );
      }
      outputFields.append( source->fields().field( groupFieldIndex ) );
    }
    else
    {
      outputFields.append( QgsField( u"group"_s, QMetaType::QString ) );
    }
  }
  outputFields.append( QgsField( "begin", orderFieldType ) );
  outputFields.append( QgsField( "end", orderFieldType ) );
 
  const bool naturalSort = parameterAsBool( parameters, u"NATURAL_SORT"_s, context );
  QCollator collator;
  collator.setNumericMode( true );
 
  Qgis::WkbType wkbType = Qgis::WkbType::LineString;
  if ( QgsWkbTypes::hasM( source->wkbType() ) )
    wkbType = QgsWkbTypes::addM( wkbType );
  if ( QgsWkbTypes::hasZ( source->wkbType() ) )
    wkbType = QgsWkbTypes::addZ( wkbType );
 
  QString dest;
  std::unique_ptr<QgsFeatureSink> sink( parameterAsSink( parameters, u"OUTPUT"_s, context, dest, outputFields, wkbType, source->sourceCrs() ) );
  if ( !sink )
    throw QgsProcessingException( invalidSinkError( parameters, u"OUTPUT"_s ) );
 
  const QString textDir = parameterAsString( parameters, u"OUTPUT_TEXT_DIR"_s, context );
  if ( !textDir.isEmpty() && !QDir().mkpath( textDir ) )
  {
    throw QgsProcessingException( QObject::tr( "Failed to create the text output directory" ) );
  }
 
  QgsDistanceArea da = QgsDistanceArea();
  da.setSourceCrs( source->sourceCrs(), context.transformContext() );
  da.setEllipsoid( context.ellipsoid() );
 
  // Store the points in a hash with the group identifier as the key
  QHash<QVariant, QVector<QPair<QVariant, QgsPoint>>> allPoints;
 
  const QgsFeatureRequest request = QgsFeatureRequest().setSubsetOfAttributes( requiredFields, source->fields() );
  QgsFeatureIterator fit = source->getFeatures( request, Qgis::ProcessingFeatureSourceFlag::SkipGeometryValidityChecks );
  QgsFeature f;
  const double totalPoints = source->featureCount() > 0 ? 100.0 / source->featureCount() : 0;
  long currentPoint = 0;
  feedback->setProgressText( QObject::tr( "Loading points…" ) );
  while ( fit.nextFeature( f ) )
  {
    if ( feedback->isCanceled() )
    {
      break;
    }
    feedback->setProgress( 0.5 * currentPoint * totalPoints );
 
    if ( f.hasGeometry() )
    {
      expressionContext.setFeature( f );
      const QVariant orderValue = orderExpression.evaluate( &expressionContext );
      const QVariant groupValue = groupExpressionString.isEmpty() ? QVariant() : groupExpression.evaluate( &expressionContext );
 
      if ( !allPoints.contains( groupValue ) )
        allPoints[groupValue] = QVector<QPair<QVariant, QgsPoint>>();
      const QgsAbstractGeometry *geom = f.geometry().constGet();
      if ( QgsWkbTypes::isMultiType( geom->wkbType() ) )
      {
        const QgsMultiPoint mp( *qgsgeometry_cast<const QgsMultiPoint *>( geom ) );
        for ( auto pit = mp.const_parts_begin(); pit != mp.const_parts_end(); ++pit )
        {
          if ( const QgsPoint *point = qgsgeometry_cast<const QgsPoint *>( *pit ) )
          {
            allPoints[groupValue] << qMakePair( orderValue, *point );
          }
        }
      }
      else
      {
        if ( const QgsPoint *point = qgsgeometry_cast<const QgsPoint *>( geom ) )
        {
          allPoints[groupValue] << qMakePair( orderValue, *point );
        }
      }
    }
    ++currentPoint;
  }
 
  int pathCount = 0;
  currentPoint = 0;
  QHashIterator<QVariant, QVector<QPair<QVariant, QgsPoint>>> hit( allPoints );
  feedback->setProgressText( QObject::tr( "Creating paths…" ) );
  while ( hit.hasNext() )
  {
    hit.next();
    if ( feedback->isCanceled() )
    {
      break;
    }
    QVector<QPair<QVariant, QgsPoint>> pairs = hit.value();
 
    if ( naturalSort )
    {
      std::stable_sort( pairs.begin(), pairs.end(), [&collator]( const QPair<QVariant, QgsPoint> &pair1, const QPair<QVariant, QgsPoint> &pair2 ) {
        return collator.compare( pair1.first.toString(), pair2.first.toString() ) < 0;
      } );
    }
    else
    {
      std::stable_sort( pairs.begin(), pairs.end(), []( const QPair<QVariant, QgsPoint> &pair1, const QPair<QVariant, QgsPoint> &pair2 ) { return qgsVariantLessThan( pair1.first, pair2.first ); } );
    }
 
 
    QVector<QgsPoint> pathPoints;
    for ( auto pit = pairs.constBegin(); pit != pairs.constEnd(); ++pit )
    {
      if ( feedback->isCanceled() )
      {
        break;
      }
      feedback->setProgress( 50 + 0.5 * currentPoint * totalPoints );
      pathPoints.append( pit->second );
      ++currentPoint;
    }
    if ( pathPoints.size() < 2 )
    {
      feedback->pushInfo( QObject::tr( "Skipping path with group %1 : insufficient vertices" ).arg( hit.key().toString() ) );
      continue;
    }
    if ( closePaths && pathPoints.size() > 2 && pathPoints.constFirst() != pathPoints.constLast() )
      pathPoints.append( pathPoints.constFirst() );
 
    QgsFeature outputFeature;
    QgsAttributes attrs;
    if ( !groupExpressionString.isEmpty() )
      attrs.append( hit.key() );
    attrs.append( pairs.first().first );
    attrs.append( pairs.last().first );
    outputFeature.setGeometry( QgsGeometry::fromPolyline( pathPoints ) );
    outputFeature.setAttributes( attrs );
    if ( !sink->addFeature( outputFeature, QgsFeatureSink::FastInsert ) )
      throw QgsProcessingException( writeFeatureError( sink.get(), parameters, u"OUTPUT"_s ) );
 
    if ( !textDir.isEmpty() )
    {
      const QString filename = QDir( textDir ).filePath( hit.key().toString() + QString( ".txt" ) );
      QFile textFile( filename );
      if ( !textFile.open( QIODevice::WriteOnly | QIODevice::Truncate ) )
        throw QgsProcessingException( QObject::tr( "Cannot open file for writing " ) + filename );
 
      QTextStream out( &textFile );
      out << QString(
               "angle=Azimuth\n"
               "heading=Coordinate_System\n"
               "dist_units=Default\n"
               "startAt=%1;%2;90\n"
               "survey=Polygonal\n"
               "[data]\n"
      )
               .arg( pathPoints.at( 0 ).x() )
               .arg( pathPoints.at( 0 ).y() );
 
      for ( int i = 1; i < pathPoints.size(); ++i )
      {
        const double angle = pathPoints.at( i - 1 ).azimuth( pathPoints.at( i ) );
        double distance = 0;
        try
        {
          distance = da.measureLine( pathPoints.at( i - 1 ), pathPoints.at( i ) );
        }
        catch ( QgsCsException & )
        {
          throw QgsProcessingException( QObject::tr( "An error occurred while calculating length" ) );
        }
        out << QString( "%1;%2;90\n" ).arg( angle ).arg( distance );
      }
    }
 
    ++pathCount;
  }
 
  sink->finalize();
 
  QVariantMap outputs;
  outputs.insert( u"OUTPUT"_s, dest );
  outputs.insert( u"NUM_PATHS"_s, pathCount );
  if ( !textDir.isEmpty() )
  {
    outputs.insert( u"OUTPUT_TEXT_DIR"_s, textDir );
  }
  return outputs;
}