23 #include <QProgressDialog> 26 #include <cpl_string.h> 27 #include <gdalwarper.h> 29 #if defined(GDAL_VERSION_NUM) && GDAL_VERSION_NUM >= 1800 30 #define TO8(x) (x).toUtf8().constData() 31 #define TO8F(x) (x).toUtf8().constData() 33 #define TO8(x) (x).toLocal8Bit().constData() 34 #define TO8F(x) QFile::encodeName( x ).constData() 39 : mFormulaString( formulaString )
40 , mOutputFile( outputFile )
41 , mOutputFormat( outputFormat )
42 , mOutputRectangle( outputExtent )
43 , mNumOutputColumns( nOutputColumns )
44 , mNumOutputRows( nOutputRows )
45 , mRasterEntries( rasterEntries )
48 mOutputCrs = mRasterEntries.
at( 0 ).raster->crs();
53 : mFormulaString( formulaString )
54 , mOutputFile( outputFile )
55 , mOutputFormat( outputFormat )
56 , mOutputRectangle( outputExtent )
57 , mOutputCrs( outputCrs )
58 , mNumOutputColumns( nOutputColumns )
59 , mNumOutputRows( nOutputRows )
60 , mRasterEntries( rasterEntries )
77 for ( ; it != mRasterEntries.
constEnd(); ++it )
82 qDeleteAll( inputBlocks );
88 if ( it->raster->crs() != mOutputCrs )
91 proj.
setCRS( it->raster->crs(), mOutputCrs );
92 proj.
setInput( it->raster->dataProvider() );
95 block = proj.
block( it->bandNumber, mOutputRectangle, mNumOutputColumns, mNumOutputRows );
99 block = it->raster->dataProvider()->block( it->bandNumber, mOutputRectangle, mNumOutputColumns, mNumOutputRows );
105 qDeleteAll( inputBlocks );
108 inputBlocks.
insert( it->ref, block );
112 GDALDriverH outputDriver = openOutputDriver();
118 GDALDatasetH outputDataset = openOutputFile( outputDriver );
120 GDALRasterBandH outputRasterBand = GDALGetRasterBand( outputDataset, 1 );
122 float outputNodataValue = -FLT_MAX;
123 GDALSetRasterNoDataValue( outputRasterBand, outputNodataValue );
134 for (
int i = 0; i < mNumOutputRows; ++i )
146 if ( calcNode->
calculate( inputBlocks, resultMatrix, i ) )
148 bool resultIsNumber = resultMatrix.
isNumber();
149 float* calcData =
new float[mNumOutputColumns];
151 for (
int j = 0; j < mNumOutputColumns; ++j )
153 calcData[j] = ( float )( resultIsNumber ? resultMatrix.
number() : resultMatrix.
data()[j] );
157 if ( GDALRasterIO( outputRasterBand, GF_Write, 0, i, mNumOutputColumns, 1, calcData, mNumOutputColumns, 1, GDT_Float32, 0, 0 ) != CE_None )
159 qWarning(
"RasterIO error!" );
174 qDeleteAll( inputBlocks );
180 GDALDeleteDataset( outputDriver,
TO8F( mOutputFile ) );
183 GDALClose( outputDataset );
185 return static_cast< int >(
Success );
189 : mNumOutputColumns( 0 )
190 , mNumOutputRows( 0 )
194 GDALDriverH QgsRasterCalculator::openOutputDriver()
196 char **driverMetadata;
199 GDALDriverH outputDriver = GDALGetDriverByName( mOutputFormat.
toLocal8Bit().
data() );
206 driverMetadata = GDALGetMetadata( outputDriver,
nullptr );
207 if ( !CSLFetchBoolean( driverMetadata, GDAL_DCAP_CREATE,
false ) )
215 GDALDatasetH QgsRasterCalculator::openOutputFile( GDALDriverH outputDriver )
218 char **papszOptions =
nullptr;
219 GDALDatasetH outputDataset = GDALCreate( outputDriver,
TO8F( mOutputFile ), mNumOutputColumns, mNumOutputRows, 1, GDT_Float32, papszOptions );
220 if ( !outputDataset )
222 return outputDataset;
226 double geotransform[6];
227 outputGeoTransform( geotransform );
228 GDALSetGeoTransform( outputDataset, geotransform );
230 return outputDataset;
233 void QgsRasterCalculator::outputGeoTransform(
double* transform )
const 235 transform[0] = mOutputRectangle.
xMinimum();
236 transform[1] = mOutputRectangle.
width() / mNumOutputColumns;
238 transform[3] = mOutputRectangle.
yMaximum();
240 transform[5] = -mOutputRectangle.
height() / mNumOutputRows;
A rectangle specified with double values.
bool calculate(QMap< QString, QgsRasterBlock * > &rasterData, QgsRasterMatrix &result, int row=-1) const
Calculates result of raster calculation (might be real matrix or single number).
void setCRS(const QgsCoordinateReferenceSystem &theSrcCRS, const QgsCoordinateReferenceSystem &theDestCRS, int srcDatumTransform=-1, int destDatumTransform=-1)
set source and destination CRS
void setMaximum(int maximum)
double yMaximum() const
Get the y maximum value (top side of rectangle)
void setNodataValue(double d)
const_iterator constEnd() const
QString toWkt() const
Returns a WKT representation of this CRS.
QgsRasterCalculator(const QString &formulaString, const QString &outputFile, const QString &outputFormat, const QgsRectangle &outputExtent, int nOutputColumns, int nOutputRows, const QVector< QgsRasterCalculatorEntry > &rasterEntries)
QgsRasterCalculator constructor.
void setValue(int progress)
QgsRasterBlock * block(int bandNo, const QgsRectangle &extent, int width, int height) override
Read block of data using given extent and size.
QByteArray toLocal8Bit() const
bool isNumber() const
Returns true if matrix is 1x1 (=scalar number)
void setPrecision(Precision precision)
int processCalculation(QProgressDialog *p=nullptr)
Starts the calculation and writes new raster.
virtual bool setInput(QgsRasterInterface *input)
Set input.
const T & at(int i) const
const_iterator constBegin() const
Class for storing a coordinate reference system (CRS)
double * data()
Returns data array (but not ownership)
static QgsRasterCalcNode * parseRasterCalcString(const QString &str, QString &parserErrorMsg)
iterator insert(const Key &key, const T &value)
double width() const
Width of the rectangle.
double xMinimum() const
Get the x minimum value (left side of rectangle)
double height() const
Height of the rectangle.
bool isEmpty() const
Returns true if block is empty, i.e.