QGIS API Documentation 3.43.0-Master (e01d6d7c4c0)
qgsalgorithmrescaleraster.cpp
Go to the documentation of this file.
1/***************************************************************************
2 qgsalgorithmrescaleraster.cpp
3 ---------------------
4 begin : July 2020
5 copyright : (C) 2020 by Alexander Bruy
6 email : alexander dot bruy 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 <limits>
19#include "math.h"
21#include "qgsrasterfilewriter.h"
22
24
25QString QgsRescaleRasterAlgorithm::name() const
26{
27 return QStringLiteral( "rescaleraster" );
28}
29
30QString QgsRescaleRasterAlgorithm::displayName() const
31{
32 return QObject::tr( "Rescale raster" );
33}
34
35QStringList QgsRescaleRasterAlgorithm::tags() const
36{
37 return QObject::tr( "raster,rescale,minimum,maximum,range" ).split( ',' );
38}
39
40QString QgsRescaleRasterAlgorithm::group() const
41{
42 return QObject::tr( "Raster analysis" );
43}
44
45QString QgsRescaleRasterAlgorithm::groupId() const
46{
47 return QStringLiteral( "rasteranalysis" );
48}
49
50QString QgsRescaleRasterAlgorithm::shortHelpString() const
51{
52 return QObject::tr( "This algorithm rescales a raster layer to a new value range, while preserving the shape "
53 "(distribution) of the raster's histogram (pixel values). Input values "
54 "are mapped using a linear interpolation from the source raster's minimum "
55 "and maximum pixel values to the destination minimum and maximum pixel range.\n\n"
56 "By default the algorithm preserves the original NoData value, but there is "
57 "an option to override it." );
58}
59
60QString QgsRescaleRasterAlgorithm::shortDescription() const
61{
62 return QObject::tr( "Rescales a raster layer to a new value range, while preserving the shape "
63 "(distribution) of the raster's histogram (pixel values)." );
64}
65
66QgsRescaleRasterAlgorithm *QgsRescaleRasterAlgorithm::createInstance() const
67{
68 return new QgsRescaleRasterAlgorithm();
69}
70
71void QgsRescaleRasterAlgorithm::initAlgorithm( const QVariantMap & )
72{
73 addParameter( new QgsProcessingParameterRasterLayer( QStringLiteral( "INPUT" ), QStringLiteral( "Input raster" ) ) );
74 addParameter( new QgsProcessingParameterBand( QStringLiteral( "BAND" ), QObject::tr( "Band number" ), 1, QStringLiteral( "INPUT" ) ) );
75 addParameter( new QgsProcessingParameterNumber( QStringLiteral( "MINIMUM" ), QObject::tr( "New minimum value" ), Qgis::ProcessingNumberParameterType::Double, 0 ) );
76 addParameter( new QgsProcessingParameterNumber( QStringLiteral( "MAXIMUM" ), QObject::tr( "New maximum value" ), Qgis::ProcessingNumberParameterType::Double, 255 ) );
77 addParameter( new QgsProcessingParameterNumber( QStringLiteral( "NODATA" ), QObject::tr( "New NoData value" ), Qgis::ProcessingNumberParameterType::Double, QVariant(), true ) );
78
79 // backwards compatibility parameter
80 // TODO QGIS 4: remove parameter and related logic
81 auto createOptsParam = std::make_unique<QgsProcessingParameterString>( QStringLiteral( "CREATE_OPTIONS" ), QObject::tr( "Creation options" ), QVariant(), false, true );
82 createOptsParam->setMetadata( QVariantMap( { { QStringLiteral( "widget_wrapper" ), QVariantMap( { { QStringLiteral( "widget_type" ), QStringLiteral( "rasteroptions" ) } } ) } } ) );
83 createOptsParam->setFlags( createOptsParam->flags() | Qgis::ProcessingParameterFlag::Hidden );
84 addParameter( createOptsParam.release() );
85
86 auto creationOptsParam = std::make_unique<QgsProcessingParameterString>( QStringLiteral( "CREATION_OPTIONS" ), QObject::tr( "Creation options" ), QVariant(), false, true );
87 creationOptsParam->setMetadata( QVariantMap( { { QStringLiteral( "widget_wrapper" ), QVariantMap( { { QStringLiteral( "widget_type" ), QStringLiteral( "rasteroptions" ) } } ) } } ) );
88 creationOptsParam->setFlags( creationOptsParam->flags() | Qgis::ProcessingParameterFlag::Advanced );
89 addParameter( creationOptsParam.release() );
90
91 addParameter( new QgsProcessingParameterRasterDestination( QStringLiteral( "OUTPUT" ), QObject::tr( "Rescaled" ) ) );
92}
93
94bool QgsRescaleRasterAlgorithm::prepareAlgorithm( const QVariantMap &parameters, QgsProcessingContext &context, QgsProcessingFeedback *feedback )
95{
96 Q_UNUSED( feedback );
97
98 QgsRasterLayer *layer = parameterAsRasterLayer( parameters, QStringLiteral( "INPUT" ), context );
99 if ( !layer )
100 throw QgsProcessingException( invalidRasterError( parameters, QStringLiteral( "INPUT" ) ) );
101
102 mBand = parameterAsInt( parameters, QStringLiteral( "BAND" ), context );
103 if ( mBand < 1 || mBand > layer->bandCount() )
104 throw QgsProcessingException( QObject::tr( "Invalid band number for BAND (%1): Valid values for input raster are 1 to %2" )
105 .arg( mBand )
106 .arg( layer->bandCount() ) );
107
108 mMinimum = parameterAsDouble( parameters, QStringLiteral( "MINIMUM" ), context );
109 mMaximum = parameterAsDouble( parameters, QStringLiteral( "MAXIMUM" ), context );
110
111 mInterface.reset( layer->dataProvider()->clone() );
112
113 mCrs = layer->crs();
114 mLayerWidth = layer->width();
115 mLayerHeight = layer->height();
116 mExtent = layer->extent();
117 if ( parameters.value( QStringLiteral( "NODATA" ) ).isValid() )
118 {
119 mNoData = parameterAsDouble( parameters, QStringLiteral( "NODATA" ), context );
120 }
121 else
122 {
123 mNoData = layer->dataProvider()->sourceNoDataValue( mBand );
124 }
125
126 if ( std::isfinite( mNoData ) )
127 {
128 // Clamp nodata to float32 range, since that's the type of the raster
129 if ( mNoData < std::numeric_limits<float>::lowest() )
130 mNoData = std::numeric_limits<float>::lowest();
131 else if ( mNoData > std::numeric_limits<float>::max() )
132 mNoData = std::numeric_limits<float>::max();
133 }
134
135 mXSize = mInterface->xSize();
136 mYSize = mInterface->ySize();
137
138 return true;
139}
140
141QVariantMap QgsRescaleRasterAlgorithm::processAlgorithm( const QVariantMap &parameters, QgsProcessingContext &context, QgsProcessingFeedback *feedback )
142{
143 feedback->pushInfo( QObject::tr( "Calculating raster minimum and maximum values…" ) );
144 const QgsRasterBandStats stats = mInterface->bandStatistics( mBand, Qgis::RasterBandStatistic::Min | Qgis::RasterBandStatistic::Max, QgsRectangle(), 0 );
145
146 feedback->pushInfo( QObject::tr( "Rescaling values…" ) );
147
148 QString creationOptions = parameterAsString( parameters, QStringLiteral( "CREATION_OPTIONS" ), context ).trimmed();
149 // handle backwards compatibility parameter CREATE_OPTIONS
150 const QString optionsString = parameterAsString( parameters, QStringLiteral( "CREATE_OPTIONS" ), context );
151 if ( !optionsString.isEmpty() )
152 creationOptions = optionsString;
153
154 const QString outputFile = parameterAsOutputLayer( parameters, QStringLiteral( "OUTPUT" ), context );
155 const QFileInfo fi( outputFile );
156 const QString outputFormat = QgsRasterFileWriter::driverForExtension( fi.suffix() );
157 auto writer = std::make_unique<QgsRasterFileWriter>( outputFile );
158 writer->setOutputProviderKey( QStringLiteral( "gdal" ) );
159 if ( !creationOptions.isEmpty() )
160 {
161 writer->setCreationOptions( creationOptions.split( '|' ) );
162 }
163
164 writer->setOutputFormat( outputFormat );
165 std::unique_ptr<QgsRasterDataProvider> provider( writer->createOneBandRaster( Qgis::DataType::Float32, mXSize, mYSize, mExtent, mCrs ) );
166 if ( !provider )
167 throw QgsProcessingException( QObject::tr( "Could not create raster output: %1" ).arg( outputFile ) );
168 if ( !provider->isValid() )
169 throw QgsProcessingException( QObject::tr( "Could not create raster output %1: %2" ).arg( outputFile, provider->error().message( QgsErrorMessage::Text ) ) );
170
171 QgsRasterDataProvider *destProvider = provider.get();
172 destProvider->setEditable( true );
173 destProvider->setNoDataValue( 1, mNoData );
174
177 const int numBlocksX = static_cast<int>( std::ceil( 1.0 * mLayerWidth / blockWidth ) );
178 const int numBlocksY = static_cast<int>( std::ceil( 1.0 * mLayerHeight / blockHeight ) );
179 const int numBlocks = numBlocksX * numBlocksY;
180
181 QgsRasterIterator iter( mInterface.get() );
182 iter.startRasterRead( mBand, mLayerWidth, mLayerHeight, mExtent );
183 int iterLeft = 0;
184 int iterTop = 0;
185 int iterCols = 0;
186 int iterRows = 0;
187 std::unique_ptr<QgsRasterBlock> inputBlock;
188 while ( iter.readNextRasterPart( mBand, iterCols, iterRows, inputBlock, iterLeft, iterTop ) )
189 {
190 auto outputBlock = std::make_unique<QgsRasterBlock>( destProvider->dataType( 1 ), iterCols, iterRows );
191 feedback->setProgress( 100 * ( ( iterTop / blockHeight * numBlocksX ) + iterLeft / blockWidth ) / numBlocks );
192
193 for ( int row = 0; row < iterRows; row++ )
194 {
195 if ( feedback->isCanceled() )
196 break;
197
198 for ( int col = 0; col < iterCols; col++ )
199 {
200 bool isNoData = false;
201 const double val = inputBlock->valueAndNoData( row, col, isNoData );
202 if ( isNoData )
203 {
204 outputBlock->setValue( row, col, mNoData );
205 }
206 else
207 {
208 const double newValue = ( ( val - stats.minimumValue ) * ( mMaximum - mMinimum ) / ( stats.maximumValue - stats.minimumValue ) ) + mMinimum;
209 outputBlock->setValue( row, col, newValue );
210 }
211 }
212 }
213 if ( !destProvider->writeBlock( outputBlock.get(), mBand, iterLeft, iterTop ) )
214 {
215 throw QgsProcessingException( QObject::tr( "Could not write raster block: %1" ).arg( destProvider->error().summary() ) );
216 }
217 }
218 destProvider->setEditable( false );
219
220 QVariantMap outputs;
221 outputs.insert( QStringLiteral( "OUTPUT" ), outputFile );
222 return outputs;
223}
224
@ Float32
Thirty two bit floating point (float)
@ Hidden
Parameter is hidden and should not be shown to users.
@ Advanced
Parameter is an advanced parameter which should be hidden from users by default.
virtual QgsError error() const
Gets current status error.
QString summary() const
Short error description, usually the first error in chain, the real error.
Definition qgserror.cpp:129
bool isCanceled() const
Tells whether the operation has been canceled already.
Definition qgsfeedback.h:53
void setProgress(double progress)
Sets the current progress for the feedback object.
Definition qgsfeedback.h:61
virtual QgsRectangle extent() const
Returns the extent of the layer.
QgsCoordinateReferenceSystem crs
Definition qgsmaplayer.h:84
Contains information about the context in which a processing algorithm is executed.
Custom exception class for processing related exceptions.
Base class for providing feedback from a processing algorithm.
virtual void pushInfo(const QString &info)
Pushes a general informational message from the algorithm.
A raster band parameter for Processing algorithms.
A numeric parameter for processing algorithms.
A raster layer destination parameter, for specifying the destination path for a raster layer created ...
A raster layer parameter for processing algorithms.
The RasterBandStats struct is a container for statistics about a single raster band.
double minimumValue
The minimum cell value in the raster band.
double maximumValue
The maximum cell value in the raster band.
Base class for raster data providers.
QgsRasterDataProvider * clone() const override=0
Clone itself, create deep copy.
virtual bool setNoDataValue(int bandNo, double noDataValue)
Set no data value on created dataset.
virtual double sourceNoDataValue(int bandNo) const
Value representing no data value.
bool writeBlock(QgsRasterBlock *block, int band, int xOffset=0, int yOffset=0)
Writes pixel data from a raster block into the provider data source.
Qgis::DataType dataType(int bandNo) const override=0
Returns data type for the band specified by number.
virtual bool setEditable(bool enabled)
Turns on/off editing mode of the provider.
static QString driverForExtension(const QString &extension)
Returns the GDAL driver name for a specified file extension.
Iterator for sequentially processing raster cells.
static const int DEFAULT_MAXIMUM_TILE_WIDTH
Default maximum tile width.
static const int DEFAULT_MAXIMUM_TILE_HEIGHT
Default maximum tile height.
Represents a raster layer.
int height() const
Returns the height of the (unclipped) raster.
int bandCount() const
Returns the number of bands in this layer.
QgsRasterDataProvider * dataProvider() override
Returns the source data provider.
int width() const
Returns the width of the (unclipped) raster.
A rectangle specified with double values.