QGIS API Documentation  3.26.3-Buenos Aires (65e4edfdad)
qgsalgorithmroundrastervalues.cpp
Go to the documentation of this file.
1 /***************************************************************************
2  qgsalgorithmroundrastervalues.cpp
3  ---------------------
4  begin : April 2020
5  copyright : (C) 2020 by Clemens Raffler
6  email : clemens dot raffler 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 "qgsrasterfilewriter.h"
20 
22 
23 QString QgsRoundRasterValuesAlgorithm::name() const
24 {
25  return QStringLiteral( "roundrastervalues" );
26 }
27 
28 QString QgsRoundRasterValuesAlgorithm::displayName() const
29 {
30  return QObject::tr( "Round raster" );
31 }
32 
33 QStringList QgsRoundRasterValuesAlgorithm::tags() const
34 {
35  return QObject::tr( "data,cells,round,truncate" ).split( ',' );
36 }
37 
38 QString QgsRoundRasterValuesAlgorithm::group() const
39 {
40  return QObject::tr( "Raster analysis" );
41 }
42 
43 QString QgsRoundRasterValuesAlgorithm::groupId() const
44 {
45  return QStringLiteral( "rasteranalysis" );
46 }
47 
48 void QgsRoundRasterValuesAlgorithm::initAlgorithm( const QVariantMap & )
49 {
50  addParameter( new QgsProcessingParameterRasterLayer( QStringLiteral( "INPUT" ), QStringLiteral( "Input raster" ) ) );
51  addParameter( new QgsProcessingParameterBand( QStringLiteral( "BAND" ), QObject::tr( "Band number" ), 1, QStringLiteral( "INPUT" ) ) );
52  addParameter( new QgsProcessingParameterEnum( QStringLiteral( "ROUNDING_DIRECTION" ), QObject::tr( "Rounding direction" ), QStringList() << QObject::tr( "Round up" ) << QObject::tr( "Round to nearest" ) << QObject::tr( "Round down" ), false, 1 ) );
53  addParameter( new QgsProcessingParameterNumber( QStringLiteral( "DECIMAL_PLACES" ), QObject::tr( "Number of decimals places" ), QgsProcessingParameterNumber::Integer, 2 ) );
54  addParameter( new QgsProcessingParameterRasterDestination( QStringLiteral( "OUTPUT" ), QObject::tr( "Output raster" ) ) );
55  std::unique_ptr< QgsProcessingParameterDefinition > baseParameter = std::make_unique< QgsProcessingParameterNumber >( QStringLiteral( "BASE_N" ), QObject::tr( "Base n for rounding to multiples of n" ), QgsProcessingParameterNumber::Integer, 10, true, 1 );
56  baseParameter->setFlags( QgsProcessingParameterDefinition::FlagAdvanced );
57  addParameter( baseParameter.release() );
58 }
59 
60 QString QgsRoundRasterValuesAlgorithm::shortHelpString() const
61 {
62  return QObject::tr( "This algorithm rounds the cell values of a raster dataset according to the specified number of decimals.\n "
63  "Alternatively, a negative number of decimal places may be used to round values to powers of a base n "
64  "(specified in the advanced parameter Base n). For example, with a Base value n of 10 and Decimal places of -1 "
65  "the algorithm rounds cell values to multiples of 10, -2 rounds to multiples of 100, and so on. Arbitrary base values "
66  "may be chosen, the algorithm applies the same multiplicative principle. Rounding cell values to multiples of "
67  "a base n may be used to generalize raster layers.\n"
68  "The algorithm preserves the data type of the input raster. Therefore byte/integer rasters can only be rounded "
69  "to multiples of a base n, otherwise a warning is raised and the raster gets copied as byte/integer raster" );
70 }
71 
72 QgsRoundRasterValuesAlgorithm *QgsRoundRasterValuesAlgorithm::createInstance() const
73 {
74  return new QgsRoundRasterValuesAlgorithm();
75 }
76 
77 bool QgsRoundRasterValuesAlgorithm::prepareAlgorithm( const QVariantMap &parameters, QgsProcessingContext &context, QgsProcessingFeedback *feedback )
78 {
79  Q_UNUSED( feedback );
80  QgsRasterLayer *inputRaster = parameterAsRasterLayer( parameters, QStringLiteral( "INPUT" ), context );
81  mDecimalPrecision = parameterAsInt( parameters, QStringLiteral( "DECIMAL_PLACES" ), context );
82  mBaseN = parameterAsInt( parameters, QStringLiteral( "BASE_N" ), context );
83  mMultipleOfBaseN = pow( mBaseN, abs( mDecimalPrecision ) );
84  mScaleFactor = std::pow( 10.0, mDecimalPrecision );
85 
86  if ( !inputRaster )
87  throw QgsProcessingException( invalidRasterError( parameters, QStringLiteral( "INPUT" ) ) );
88 
89  mBand = parameterAsInt( parameters, QStringLiteral( "BAND" ), context );
90  if ( mBand < 1 || mBand > inputRaster->bandCount() )
91  throw QgsProcessingException( QObject::tr( "Invalid band number for BAND (%1): Valid values for input raster are 1 to %2" ).arg( mBand ).arg( inputRaster->bandCount() ) );
92 
93  mRoundingDirection = parameterAsEnum( parameters, QStringLiteral( "ROUNDING_DIRECTION" ), context );
94 
95  mInterface.reset( inputRaster->dataProvider()->clone() );
96  mDataType = mInterface->dataType( mBand );
97 
98  switch ( mDataType )
99  {
105  mIsInteger = true;
106  if ( mDecimalPrecision > -1 )
107  feedback->reportError( QObject::tr( "Input raster is of byte or integer type. The cell values cannot be rounded and will be output using the same data type." ), false );
108  break;
109  default:
110  mIsInteger = false;
111  break;
112  }
113 
114  mInputNoDataValue = inputRaster->dataProvider()->sourceNoDataValue( mBand );
115  mExtent = inputRaster->extent();
116  mLayerWidth = inputRaster->width();
117  mLayerHeight = inputRaster->height();
118  mCrs = inputRaster->crs();
119  mNbCellsXProvider = mInterface->xSize();
120  mNbCellsYProvider = mInterface->ySize();
121  return true;
122 }
123 
124 QVariantMap QgsRoundRasterValuesAlgorithm::processAlgorithm( const QVariantMap &parameters, QgsProcessingContext &context, QgsProcessingFeedback *feedback )
125 {
126  //prepare output dataset
127  const QString outputFile = parameterAsOutputLayer( parameters, QStringLiteral( "OUTPUT" ), context );
128  const QFileInfo fi( outputFile );
129  const QString outputFormat = QgsRasterFileWriter::driverForExtension( fi.suffix() );
130  std::unique_ptr< QgsRasterFileWriter > writer = std::make_unique< QgsRasterFileWriter >( outputFile );
131  writer->setOutputProviderKey( QStringLiteral( "gdal" ) );
132  writer->setOutputFormat( outputFormat );
133  std::unique_ptr<QgsRasterDataProvider > provider( writer->createOneBandRaster( mInterface->dataType( mBand ), mNbCellsXProvider, mNbCellsYProvider, mExtent, mCrs ) );
134  if ( !provider )
135  throw QgsProcessingException( QObject::tr( "Could not create raster output: %1" ).arg( outputFile ) );
136  if ( !provider->isValid() )
137  throw QgsProcessingException( QObject::tr( "Could not create raster output %1: %2" ).arg( outputFile, provider->error().message( QgsErrorMessage::Text ) ) );
138 
139  //prepare output provider
140  QgsRasterDataProvider *destinationRasterProvider;
141  destinationRasterProvider = provider.get();
142  destinationRasterProvider->setEditable( true );
143  destinationRasterProvider->setNoDataValue( 1, mInputNoDataValue );
144 
147  const int nbBlocksWidth = static_cast< int >( std::ceil( 1.0 * mLayerWidth / maxWidth ) );
148  const int nbBlocksHeight = static_cast< int >( std::ceil( 1.0 * mLayerHeight / maxHeight ) );
149  const int nbBlocks = nbBlocksWidth * nbBlocksHeight;
150 
151  QgsRasterIterator iter( mInterface.get() );
152  iter.startRasterRead( mBand, mLayerWidth, mLayerHeight, mExtent );
153  int iterLeft = 0;
154  int iterTop = 0;
155  int iterCols = 0;
156  int iterRows = 0;
157  std::unique_ptr< QgsRasterBlock > analysisRasterBlock;
158  while ( iter.readNextRasterPart( mBand, iterCols, iterRows, analysisRasterBlock, iterLeft, iterTop ) )
159  {
160  if ( feedback )
161  feedback->setProgress( 100 * ( ( iterTop / maxHeight * nbBlocksWidth ) + iterLeft / maxWidth ) / nbBlocks );
162  if ( mIsInteger && mDecimalPrecision > -1 )
163  {
164  //nothing to round, just write raster block
165  analysisRasterBlock->setNoDataValue( mInputNoDataValue );
166  destinationRasterProvider->writeBlock( analysisRasterBlock.get(), mBand, iterLeft, iterTop );
167  }
168  else
169  {
170  for ( int row = 0; row < iterRows; row++ )
171  {
172  if ( feedback && feedback->isCanceled() )
173  break;
174  for ( int column = 0; column < iterCols; column++ )
175  {
176  bool isNoData = false;
177  const double val = analysisRasterBlock->valueAndNoData( row, column, isNoData );
178  if ( isNoData )
179  {
180  analysisRasterBlock->setValue( row, column, mInputNoDataValue );
181  }
182  else
183  {
184  double roundedVal = mInputNoDataValue;
185  if ( mRoundingDirection == 0 && mDecimalPrecision < 0 )
186  {
187  roundedVal = roundUpBaseN( val );
188  }
189  else if ( mRoundingDirection == 0 && mDecimalPrecision > -1 )
190  {
191  const double m = ( val < 0.0 ) ? -1.0 : 1.0;
192  roundedVal = roundUp( val, m );
193  }
194  else if ( mRoundingDirection == 1 && mDecimalPrecision < 0 )
195  {
196  roundedVal = roundNearestBaseN( val );
197  }
198  else if ( mRoundingDirection == 1 && mDecimalPrecision > -1 )
199  {
200  const double m = ( val < 0.0 ) ? -1.0 : 1.0;
201  roundedVal = roundNearest( val, m );
202  }
203  else if ( mRoundingDirection == 2 && mDecimalPrecision < 0 )
204  {
205  roundedVal = roundDownBaseN( val );
206  }
207  else
208  {
209  const double m = ( val < 0.0 ) ? -1.0 : 1.0;
210  roundedVal = roundDown( val, m );
211  }
212  //integer values get automatically cast to double when reading and back to int when writing
213  analysisRasterBlock->setValue( row, column, roundedVal );
214  }
215  }
216  }
217  destinationRasterProvider->writeBlock( analysisRasterBlock.get(), mBand, iterLeft, iterTop );
218  }
219  }
220  destinationRasterProvider->setEditable( false );
221 
222  QVariantMap outputs;
223  outputs.insert( QStringLiteral( "OUTPUT" ), outputFile );
224  return outputs;
225 }
226 
227 double QgsRoundRasterValuesAlgorithm::roundNearest( double value, double m )
228 {
229  return ( std::round( value * m * mScaleFactor ) / mScaleFactor ) * m;
230 }
231 
232 double QgsRoundRasterValuesAlgorithm::roundUp( double value, double m )
233 {
234  return ( std::ceil( value * m * mScaleFactor ) / mScaleFactor ) * m;
235 }
236 
237 double QgsRoundRasterValuesAlgorithm::roundDown( double value, double m )
238 {
239  return ( std::floor( value * m * mScaleFactor ) / mScaleFactor ) * m;
240 }
241 
242 
243 double QgsRoundRasterValuesAlgorithm::roundNearestBaseN( double value )
244 {
245  return static_cast<double>( mMultipleOfBaseN * round( value / mMultipleOfBaseN ) );
246 }
247 
248 double QgsRoundRasterValuesAlgorithm::roundUpBaseN( double value )
249 {
250  return static_cast<double>( mMultipleOfBaseN * ceil( value / mMultipleOfBaseN ) );
251 }
252 
253 double QgsRoundRasterValuesAlgorithm::roundDownBaseN( double value )
254 {
255  return static_cast<double>( mMultipleOfBaseN * floor( value / mMultipleOfBaseN ) );
256 }
257 
258 
259 
260 
QgsMapLayer::crs
QgsCoordinateReferenceSystem crs
Definition: qgsmaplayer.h:79
QgsFeedback::setProgress
void setProgress(double progress)
Sets the current progress for the feedback object.
Definition: qgsfeedback.h:76
QgsRasterLayer::bandCount
int bandCount() const
Returns the number of bands in this layer.
Definition: qgsrasterlayer.cpp:240
Qgis::DataType::UInt32
@ UInt32
Thirty two bit unsigned integer (quint32)
QgsProcessingParameterNumber
A numeric parameter for processing algorithms.
Definition: qgsprocessingparameters.h:2179
QgsProcessingFeedback
Base class for providing feedback from a processing algorithm.
Definition: qgsprocessingfeedback.h:37
QgsRasterFileWriter::driverForExtension
static QString driverForExtension(const QString &extension)
Returns the GDAL driver name for a specified file extension.
Definition: qgsrasterfilewriter.cpp:1066
QgsProcessingFeedback::reportError
virtual void reportError(const QString &error, bool fatalError=false)
Reports that the algorithm encountered an error while executing.
Definition: qgsprocessingfeedback.cpp:59
QgsProcessingParameterRasterDestination
A raster layer destination parameter, for specifying the destination path for a raster layer created ...
Definition: qgsprocessingparameters.h:3390
Qgis::DataType::Byte
@ Byte
Eight bit unsigned integer (quint8)
QgsFeedback::isCanceled
bool isCanceled() const SIP_HOLDGIL
Tells whether the operation has been canceled already.
Definition: qgsfeedback.h:67
QgsProcessingParameterDefinition::FlagAdvanced
@ FlagAdvanced
Parameter is an advanced parameter which should be hidden from users by default.
Definition: qgsprocessingparameters.h:451
Qgis::DataType::Int32
@ Int32
Thirty two bit signed integer (qint32)
QgsRasterDataProvider::setEditable
virtual bool setEditable(bool enabled)
Turns on/off editing mode of the provider.
Definition: qgsrasterdataprovider.h:479
qgsrasterfilewriter.h
QgsRasterDataProvider::clone
QgsRasterDataProvider * clone() const override=0
Clone itself, create deep copy.
QgsRasterLayer::width
int width() const
Returns the width of the (unclipped) raster.
Definition: qgsrasterlayer.cpp:2618
QgsProcessingContext
Contains information about the context in which a processing algorithm is executed.
Definition: qgsprocessingcontext.h:46
QgsRasterDataProvider::writeBlock
bool writeBlock(QgsRasterBlock *block, int band, int xOffset=0, int yOffset=0)
Writes pixel data from a raster block into the provider data source.
Definition: qgsrasterdataprovider.cpp:357
QgsRasterLayer::height
int height() const
Returns the height of the (unclipped) raster.
Definition: qgsrasterlayer.cpp:2624
QgsMapLayer::extent
virtual QgsRectangle extent() const
Returns the extent of the layer.
Definition: qgsmaplayer.cpp:305
QgsRasterIterator::DEFAULT_MAXIMUM_TILE_HEIGHT
static const int DEFAULT_MAXIMUM_TILE_HEIGHT
Default maximum tile height.
Definition: qgsrasteriterator.h:170
Qgis::DataType::Int16
@ Int16
Sixteen bit signed integer (qint16)
QgsProcessingParameterRasterLayer
A raster layer parameter for processing algorithms.
Definition: qgsprocessingparameters.h:2495
QgsRasterDataProvider::sourceNoDataValue
virtual double sourceNoDataValue(int bandNo) const
Value representing no data value.
Definition: qgsrasterdataprovider.h:250
QgsRasterLayer
Represents a raster layer.
Definition: qgsrasterlayer.h:76
QgsRasterIterator
Iterator for sequentially processing raster cells.
Definition: qgsrasteriterator.h:34
Qgis::DataType::UInt16
@ UInt16
Sixteen bit unsigned integer (quint16)
QgsProcessingParameterNumber::Integer
@ Integer
Integer values.
Definition: qgsprocessingparameters.h:2186
QgsRasterDataProvider::setNoDataValue
virtual bool setNoDataValue(int bandNo, double noDataValue)
Set no data value on created dataset.
Definition: qgsrasterdataprovider.h:527
QgsRasterIterator::DEFAULT_MAXIMUM_TILE_WIDTH
static const int DEFAULT_MAXIMUM_TILE_WIDTH
Default maximum tile width.
Definition: qgsrasteriterator.h:167
QgsErrorMessage::Text
@ Text
Definition: qgserror.h:38
QgsProcessingParameterBand
A raster band parameter for Processing algorithms.
Definition: qgsprocessingparameters.h:3548
qgsalgorithmroundrastervalues.h
QgsProcessingParameterEnum
An enum based parameter for processing algorithms, allowing for selection from predefined values.
Definition: qgsprocessingparameters.h:2540
QgsProcessingException
Custom exception class for processing related exceptions.
Definition: qgsexception.h:82
QgsRasterDataProvider
Base class for raster data providers.
Definition: qgsrasterdataprovider.h:88
QgsRasterLayer::dataProvider
QgsRasterDataProvider * dataProvider() override
Returns the source data provider.
Definition: qgsrasterlayer.cpp:257