00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018 #include "qgsninecellfilter.h"
00019 #include "cpl_string.h"
00020 #include <QProgressDialog>
00021
00022 #if defined(GDAL_VERSION_NUM) && GDAL_VERSION_NUM >= 1800
00023 #define TO8(x) (x).toUtf8().constData()
00024 #else
00025 #define TO8(x) (x).toLocal8Bit().constData()
00026 #endif
00027
00028 QgsNineCellFilter::QgsNineCellFilter( const QString& inputFile, const QString& outputFile, const QString& outputFormat ): \
00029 mInputFile( inputFile ), mOutputFile( outputFile ), mOutputFormat( outputFormat ), mCellSizeX( -1 ), mCellSizeY( -1 ), mInputNodataValue( -1 ), mOutputNodataValue( -1 )
00030 {
00031
00032 }
00033
00034 QgsNineCellFilter::QgsNineCellFilter()
00035 {
00036
00037 }
00038
00039 QgsNineCellFilter::~QgsNineCellFilter()
00040 {
00041
00042 }
00043
00044 int QgsNineCellFilter::processRaster( QProgressDialog* p )
00045 {
00046 GDALAllRegister();
00047
00048
00049 int xSize, ySize;
00050 GDALDatasetH inputDataset = openInputFile( xSize, ySize );
00051 if ( inputDataset == NULL )
00052 {
00053 return 1;
00054 }
00055
00056
00057 GDALDriverH outputDriver = openOutputDriver();
00058 if ( outputDriver == 0 )
00059 {
00060 return 2;
00061 }
00062
00063 GDALDatasetH outputDataset = openOutputFile( inputDataset, outputDriver );
00064 if ( outputDataset == NULL )
00065 {
00066 return 3;
00067 }
00068
00069
00070 GDALRasterBandH rasterBand = GDALGetRasterBand( inputDataset, 1 );
00071 if ( rasterBand == NULL )
00072 {
00073 GDALClose( inputDataset );
00074 GDALClose( outputDataset );
00075 return 4;
00076 }
00077 mInputNodataValue = GDALGetRasterNoDataValue( rasterBand, NULL );
00078
00079 GDALRasterBandH outputRasterBand = GDALGetRasterBand( outputDataset, 1 );
00080 if ( outputRasterBand == NULL )
00081 {
00082 GDALClose( inputDataset );
00083 GDALClose( outputDataset );
00084 return 5;
00085 }
00086
00087 GDALSetRasterNoDataValue( outputRasterBand, -9999 );
00088 mOutputNodataValue = GDALGetRasterNoDataValue( outputRasterBand, NULL );
00089
00090 if ( ySize < 3 )
00091 {
00092 GDALClose( inputDataset );
00093 GDALClose( outputDataset );
00094 return 6;
00095 }
00096
00097
00098 float* scanLine1 = ( float * ) CPLMalloc( sizeof( float ) * xSize );
00099 float* scanLine2 = ( float * ) CPLMalloc( sizeof( float ) * xSize );
00100 float* scanLine3 = ( float * ) CPLMalloc( sizeof( float ) * xSize );
00101
00102 float* resultLine = ( float * ) CPLMalloc( sizeof( float ) * xSize );
00103
00104 if ( p )
00105 {
00106 p->setMaximum( ySize );
00107 }
00108
00109
00110 for ( int i = 0; i < ySize; ++i )
00111 {
00112 if ( p )
00113 {
00114 p->setValue( i );
00115 }
00116
00117 if ( p && p->wasCanceled() )
00118 {
00119 break;
00120 }
00121
00122 if ( i == 0 )
00123 {
00124
00125 for ( int a = 0; a < xSize; ++a )
00126 {
00127 scanLine1[a] = mInputNodataValue;
00128 }
00129 GDALRasterIO( rasterBand, GF_Read, 0, 0, xSize, 1, scanLine2, xSize, 1, GDT_Float32, 0, 0 );
00130 }
00131 else
00132 {
00133
00134 CPLFree( scanLine1 );
00135 scanLine1 = scanLine2;
00136 scanLine2 = scanLine3;
00137 scanLine3 = ( float * ) CPLMalloc( sizeof( float ) * xSize );
00138 }
00139
00140 if ( i == ySize - 1 )
00141 {
00142 for ( int a = 0; a < xSize; ++a )
00143 {
00144 scanLine3[a] = mInputNodataValue;
00145 }
00146 }
00147 else
00148 {
00149 GDALRasterIO( rasterBand, GF_Read, 0, i + 1, xSize, 1, scanLine3, xSize, 1, GDT_Float32, 0, 0 );
00150 }
00151
00152 for ( int j = 0; j < xSize; ++j )
00153 {
00154 if ( j == 0 )
00155 {
00156 resultLine[j] = processNineCellWindow( &mInputNodataValue, &scanLine1[j], &scanLine1[j+1], &mInputNodataValue, &scanLine2[j], \
00157 &scanLine2[j+1], &mInputNodataValue, &scanLine3[j], &scanLine3[j+1] );
00158 }
00159 else if ( j == xSize - 1 )
00160 {
00161 resultLine[j] = processNineCellWindow( &scanLine1[j-1], &scanLine1[j], &mInputNodataValue, &scanLine2[j-1], &scanLine2[j], \
00162 &mInputNodataValue, &scanLine3[j-1], &scanLine3[j], &mInputNodataValue );
00163 }
00164 else
00165 {
00166 resultLine[j] = processNineCellWindow( &scanLine1[j-1], &scanLine1[j], &scanLine1[j+1], &scanLine2[j-1], &scanLine2[j], \
00167 &scanLine2[j+1], &scanLine3[j-1], &scanLine3[j], &scanLine3[j+1] );
00168 }
00169 }
00170
00171 GDALRasterIO( outputRasterBand, GF_Write, 0, i, xSize, 1, resultLine, xSize, 1, GDT_Float32, 0, 0 );
00172 }
00173
00174 if ( p )
00175 {
00176 p->setValue( ySize );
00177 }
00178
00179 CPLFree( resultLine );
00180 CPLFree( scanLine1 );
00181 CPLFree( scanLine2 );
00182 CPLFree( scanLine3 );
00183
00184 GDALClose( inputDataset );
00185
00186 if ( p && p->wasCanceled() )
00187 {
00188
00189 GDALDeleteDataset( outputDriver, mOutputFile.toLocal8Bit().data() );
00190 return 7;
00191 }
00192 GDALClose( outputDataset );
00193
00194 return 0;
00195 }
00196
00197 GDALDatasetH QgsNineCellFilter::openInputFile( int& nCellsX, int& nCellsY )
00198 {
00199 GDALDatasetH inputDataset = GDALOpen( TO8( mInputFile ), GA_ReadOnly );
00200 if ( inputDataset != NULL )
00201 {
00202 nCellsX = GDALGetRasterXSize( inputDataset );
00203 nCellsY = GDALGetRasterYSize( inputDataset );
00204
00205
00206 if ( GDALGetRasterCount( inputDataset ) < 1 )
00207 {
00208 GDALClose( inputDataset );
00209 return NULL;
00210 }
00211 }
00212 return inputDataset;
00213 }
00214
00215 GDALDriverH QgsNineCellFilter::openOutputDriver()
00216 {
00217 char **driverMetadata;
00218
00219
00220 GDALDriverH outputDriver = GDALGetDriverByName( mOutputFormat.toLocal8Bit().data() );
00221
00222 if ( outputDriver == NULL )
00223 {
00224 return outputDriver;
00225 }
00226
00227 driverMetadata = GDALGetMetadata( outputDriver, NULL );
00228 if ( !CSLFetchBoolean( driverMetadata, GDAL_DCAP_CREATE, false ) )
00229 {
00230 return NULL;
00231 }
00232
00233 return outputDriver;
00234 }
00235
00236 GDALDatasetH QgsNineCellFilter::openOutputFile( GDALDatasetH inputDataset, GDALDriverH outputDriver )
00237 {
00238 if ( inputDataset == NULL )
00239 {
00240 return NULL;
00241 }
00242
00243 int xSize = GDALGetRasterXSize( inputDataset );
00244 int ySize = GDALGetRasterYSize( inputDataset );;
00245
00246
00247 char **papszOptions = NULL;
00248 GDALDatasetH outputDataset = GDALCreate( outputDriver, mOutputFile.toLocal8Bit().data(), xSize, ySize, 1, GDT_Float32, papszOptions );
00249 if ( outputDataset == NULL )
00250 {
00251 return outputDataset;
00252 }
00253
00254
00255 double geotransform[6];
00256 if ( GDALGetGeoTransform( inputDataset, geotransform ) != CE_None )
00257 {
00258 GDALClose( outputDataset );
00259 return NULL;
00260 }
00261 GDALSetGeoTransform( outputDataset, geotransform );
00262
00263
00264 mCellSizeX = geotransform[1];
00265 if ( mCellSizeX < 0 )
00266 {
00267 mCellSizeX = -mCellSizeX;
00268 }
00269 mCellSizeY = geotransform[5];
00270 if ( mCellSizeY < 0 )
00271 {
00272 mCellSizeY = -mCellSizeY;
00273 }
00274
00275 const char* projection = GDALGetProjectionRef( inputDataset );
00276 GDALSetProjection( outputDataset, projection );
00277
00278 return outputDataset;
00279 }
00280