Skip to content
Permalink
Browse files

[FEATURE]: relief and hillshade function, modification of raster terr…

…ain structure
  • Loading branch information
mhugent committed Jan 12, 2012
1 parent 92de863 commit bd273f62c3f73b2e8d1a9f387b0e036fe437ce9b
@@ -26,9 +26,11 @@ SET(QGIS_ANALYSIS_SRCS
raster/qgsninecellfilter.cpp
raster/qgsruggednessfilter.cpp
raster/qgsderivativefilter.cpp
raster/qgshillshadefilter.cpp
raster/qgsslopefilter.cpp
raster/qgsaspectfilter.cpp
raster/qgstotalcurvaturefilter.cpp
raster/qgsrelief.cpp
raster/qgsrastercalcnode.cpp
raster/qgsrastercalculator.cpp
raster/qgsrastermatrix.cpp
@@ -27,13 +27,11 @@ class ANALYSIS_EXPORT QgsAspectFilter: public QgsDerivativeFilter
QgsAspectFilter( const QString& inputFile, const QString& outputFile, const QString& outputFormat );
~QgsAspectFilter();

protected:
protected:
/**Calculates output value from nine input values. The input values and the output value can be equal to the
/**Calculates output value from nine input values. The input values and the output value can be equal to the \
nodata value if not present or outside of the border. Must be implemented by subclasses*/
float processNineCellWindow( float* x11, float* x21, float* x31,
float* x12, float* x22, float* x32,
float* x13, float* x23, float* x33 );
float processNineCellWindow( float* x11, float* x21, float* x31, \
float* x12, float* x22, float* x32, float* x13, float* x23, float* x33 );

};

#endif // QGSASPECTFILTER_H
@@ -17,8 +17,8 @@

#include "qgsderivativefilter.h"

QgsDerivativeFilter::QgsDerivativeFilter( const QString& inputFile, const QString& outputFile, const QString& outputFormat )
: QgsNineCellFilter( inputFile, outputFile, outputFormat )
QgsDerivativeFilter::QgsDerivativeFilter( const QString& inputFile, const QString& outputFile, const QString& outputFormat ): \
QgsNineCellFilter( inputFile, outputFile, outputFormat )
{

}
@@ -92,7 +92,7 @@ float QgsDerivativeFilter::calcFirstDerX( float* x11, float* x21, float* x31, fl
return mOutputNodataValue;
}

return sum / ( weight * mCellSizeX );
return sum / ( weight * mCellSizeX * mZFactor );
}

float QgsDerivativeFilter::calcFirstDerY( float* x11, float* x21, float* x31, float* x12, float* x22, float* x32, float* x13, float* x23, float* x33 )
@@ -159,7 +159,7 @@ float QgsDerivativeFilter::calcFirstDerY( float* x11, float* x21, float* x31, fl
return mOutputNodataValue;
}

return sum / ( weight * mCellSizeY );
return sum / ( weight * mCellSizeY * mZFactor );
}


@@ -27,9 +27,8 @@ class QgsDerivativeFilter: public QgsNineCellFilter
QgsDerivativeFilter( const QString& inputFile, const QString& outputFile, const QString& outputFormat );
virtual ~QgsDerivativeFilter();
//to be implemented by subclasses
virtual float processNineCellWindow( float* x11, float* x21, float* x31,
float* x12, float* x22, float* x32,
float* x13, float* x23, float* x33 ) = 0;
virtual float processNineCellWindow( float* x11, float* x21, float* x31, float* x12, float* x22, \
float* x32, float* x13, float* x23, float* x33 ) = 0;

protected:
/**Calculates the first order derivative in x-direction according to Horn (1981)*/
@@ -0,0 +1,54 @@
/***************************************************************************
qgshillshadefilter.h - description
--------------------------------
begin : September 26th, 2011
copyright : (C) 2011 by Marco Hugentobler
email : marco dot hugentobler at sourcepole dot ch
***************************************************************************/

/***************************************************************************
* *
* 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 "qgshillshadefilter.h"

QgsHillshadeFilter::QgsHillshadeFilter( const QString& inputFile, const QString& outputFile, const QString& outputFormat, double lightAzimuth,
double lightAngle): \
QgsDerivativeFilter( inputFile, outputFile, outputFormat ), mLightAzimuth( lightAzimuth ), mLightAngle( lightAngle )
{
}

QgsHillshadeFilter::~QgsHillshadeFilter()
{
}

float QgsHillshadeFilter::processNineCellWindow( float* x11, float* x21, float* x31, \
float* x12, float* x22, float* x32, float* x13, float* x23, float* x33 )
{
float derX = calcFirstDerX( x11, x21, x31, x12, x22, x32, x13, x23, x33 );
float derY = calcFirstDerY( x11, x21, x31, x12, x22, x32, x13, x23, x33 );

if ( derX == mOutputNodataValue || derY == mOutputNodataValue )
{
return mOutputNodataValue;
}

float zenith_rad = mLightAngle * M_PI / 180.0;
float slope_rad = atan( sqrt( derX * derX + derY * derY ) );
float azimuth_rad = mLightAzimuth * M_PI / 180.0;
float aspect_rad = 0;
if( derX == 0 && derY == 0 ) //aspect undefined, take a neutral value. Better solutions?
{
aspect_rad = azimuth_rad / 2.0;
}
else
{
aspect_rad = M_PI + atan2( derX, derY );
}
return qMax( 0.0, 255.0 * ( ( cos( zenith_rad ) * cos( slope_rad ) ) + ( sin( zenith_rad ) * sin( slope_rad ) * cos( azimuth_rad - aspect_rad ) ) ) );
}
@@ -0,0 +1,45 @@
/***************************************************************************
qgshillshadefilter.h - description
--------------------------------
begin : September 26th, 2011
copyright : (C) 2011 by Marco Hugentobler
email : marco dot hugentobler at sourcepole dot ch
***************************************************************************/

/***************************************************************************
* *
* 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. *
* *
***************************************************************************/

#ifndef QGSHILLSHADEFILTER_H
#define QGSHILLSHADEFILTER_H

#include "qgsderivativefilter.h"

class ANALYSIS_EXPORT QgsHillshadeFilter: public QgsDerivativeFilter
{
public:
QgsHillshadeFilter( const QString& inputFile, const QString& outputFile, const QString& outputFormat, double lightAzimuth = 300,
double lightAngle = 40 );
~QgsHillshadeFilter();

/**Calculates output value from nine input values. The input values and the output value can be equal to the \
nodata value if not present or outside of the border. Must be implemented by subclasses*/
float processNineCellWindow( float* x11, float* x21, float* x31, \
float* x12, float* x22, float* x32, float* x13, float* x23, float* x33 );

float lightAzimuth() const { return mLightAzimuth; }
void setLightAzimuth( float azimuth ){ mLightAzimuth = azimuth; }
float lightAngle() const { return mLightAngle; }
void setLightAngle( float angle ){ mLightAngle = angle; }

private:
float mLightAzimuth;
float mLightAngle;
};

#endif // QGSHILLSHADEFILTER_H
@@ -25,8 +25,9 @@
#define TO8(x) (x).toLocal8Bit().constData()
#endif

QgsNineCellFilter::QgsNineCellFilter( const QString& inputFile, const QString& outputFile, const QString& outputFormat )
: mInputFile( inputFile ), mOutputFile( outputFile ), mOutputFormat( outputFormat ), mCellSizeX( -1 ), mCellSizeY( -1 ), mInputNodataValue( -1 ), mOutputNodataValue( -1 )
QgsNineCellFilter::QgsNineCellFilter( const QString& inputFile, const QString& outputFile, const QString& outputFormat ): \
mInputFile( inputFile ), mOutputFile( outputFile ), mOutputFormat( outputFormat ), mCellSizeX( -1 ), mCellSizeY( -1 ), \
mInputNodataValue( -1 ), mOutputNodataValue( -1 ), mZFactor( 1.0/*111120*/ )
{

}
@@ -153,17 +154,17 @@ int QgsNineCellFilter::processRaster( QProgressDialog* p )
{
if ( j == 0 )
{
resultLine[j] = processNineCellWindow( &mInputNodataValue, &scanLine1[j], &scanLine1[j+1], &mInputNodataValue, &scanLine2[j],
resultLine[j] = processNineCellWindow( &mInputNodataValue, &scanLine1[j], &scanLine1[j+1], &mInputNodataValue, &scanLine2[j], \
&scanLine2[j+1], &mInputNodataValue, &scanLine3[j], &scanLine3[j+1] );
}
else if ( j == xSize - 1 )
{
resultLine[j] = processNineCellWindow( &scanLine1[j-1], &scanLine1[j], &mInputNodataValue, &scanLine2[j-1], &scanLine2[j],
resultLine[j] = processNineCellWindow( &scanLine1[j-1], &scanLine1[j], &mInputNodataValue, &scanLine2[j-1], &scanLine2[j], \
&mInputNodataValue, &scanLine3[j-1], &scanLine3[j], &mInputNodataValue );
}
else
{
resultLine[j] = processNineCellWindow( &scanLine1[j-1], &scanLine1[j], &scanLine1[j+1], &scanLine2[j-1], &scanLine2[j],
resultLine[j] = processNineCellWindow( &scanLine1[j-1], &scanLine1[j], &scanLine1[j+1], &scanLine2[j-1], &scanLine2[j], \
&scanLine2[j+1], &scanLine3[j-1], &scanLine3[j], &scanLine3[j+1] );
}
}
@@ -23,8 +23,8 @@

class QProgressDialog;

/**Base class for raster analysis methods that work with a 3x3 cell filter and calculate the value of each cell based on
the cell value and the eight neighbour cells. Common examples are slope and aspect calculation in DEMs. Subclasses only implement
/**Base class for raster analysis methods that work with a 3x3 cell filter and calculate the value of each cell based on \
the cell value and the eight neighbour cells. Common examples are slope and aspect calculation in DEMs. Subclasses only implement \
the method that calculates the new value from the nine values. Everything else (reading file, writing file) is done by this subclass*/

class ANALYSIS_EXPORT QgsNineCellFilter
@@ -38,6 +38,24 @@ class ANALYSIS_EXPORT QgsNineCellFilter
@return 0 in case of success*/
int processRaster( QProgressDialog* p );

double cellSizeX() const { return mCellSizeX; }
void setCellSizeX( double size ) { mCellSizeX = size; }
double cellSizeY() const { return mCellSizeY; }
void setCellSizeY( double size ) { mCellSizeY = size; }

double zFactor() const { return mZFactor; }
void setZFactor( double factor ) { mZFactor = factor; }

double inputNodataValue() const { return mInputNodataValue; }
void setInputNodataValue( double value ){ mInputNodataValue = value; }
double outputNodataValue() const { return mOutputNodataValue; }
void setOutputNodataValue( double value ){ mOutputNodataValue = value; }

/**Calculates output value from nine input values. The input values and the output value can be equal to the \
nodata value if not present or outside of the border. Must be implemented by subclasses*/
virtual float processNineCellWindow( float* x11, float* x21, float* x31, \
float* x12, float* x22, float* x32, float* x13, float* x23, float* x33 ) = 0;

private:
//default constructor forbidden. We need input file, output file and format obligatory
QgsNineCellFilter();
@@ -52,11 +70,6 @@ class ANALYSIS_EXPORT QgsNineCellFilter
GDALDatasetH openOutputFile( GDALDatasetH inputDataset, GDALDriverH outputDriver );

protected:
/**Calculates output value from nine input values. The input values and the output value can be equal to the
nodata value if not present or outside of the border. Must be implemented by subclasses*/
virtual float processNineCellWindow( float* x11, float* x21, float* x31,
float* x12, float* x22, float* x32,
float* x13, float* x23, float* x33 ) = 0;

QString mInputFile;
QString mOutputFile;
@@ -68,6 +81,8 @@ class ANALYSIS_EXPORT QgsNineCellFilter
float mInputNodataValue;
/**The nodata value of the output layer*/
float mOutputNodataValue;
/**Scale factor for z-value if x-/y- units are different to z-units (111120 for degree->meters and 370400 for degree->feet)*/
double mZFactor;
};

#endif // QGSNINECELLFILTER_H

0 comments on commit bd273f6

Please sign in to comment.
You can’t perform that action at this time.