qgsderivativefilter.cpp

Go to the documentation of this file.

/***************************************************************************
                          qgsderivativefilter.cpp  -  description
                          -----------------------
    begin                : August 7th, 2009
    copyright            : (C) 2009 by Marco Hugentobler
    email                : marco dot hugentobler at karto dot baug dot ethz 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 "qgsderivativefilter.h"
 
QgsDerivativeFilter::QgsDerivativeFilter( const QString &inputFile, const QString &outputFile, const QString &outputFormat )
  : QgsNineCellFilter( inputFile, outputFile, outputFormat )
{
}
 
float QgsDerivativeFilter::calcFirstDerX( float x11, float x21, float x31, float x12, float x22, float x32, float x13, float x23, float x33 ) const
{
  //the basic formula would be simple, but we need to test for nodata values...
  //return (( (*x31 - *x11) + 2 * (*x32 - *x12) + (*x33 - *x13) ) / (8 * mCellSizeX));
 
  int weight = 0;
  double sum = 0;
 
  //first row
  if ( x31 != mInputNodataValue && x11 != mInputNodataValue ) //the normal case
  {
    sum += ( x31 - x11 );
    weight += 2;
  }
  else if ( x31 == mInputNodataValue && x11 != mInputNodataValue && x21 != mInputNodataValue ) //probably 3x3 window is at the border
  {
    sum += ( x21 - x11 );
    weight += 1;
  }
  else if ( x11 == mInputNodataValue && x31 != mInputNodataValue && x21 != mInputNodataValue ) //probably 3x3 window is at the border
  {
    sum += ( x31 - x21 );
    weight += 1;
  }
 
  //second row
  if ( x32 != mInputNodataValue && x12 != mInputNodataValue ) //the normal case
  {
    sum += 2 * ( x32 - x12 );
    weight += 4;
  }
  else if ( x32 == mInputNodataValue && x12 != mInputNodataValue && x22 != mInputNodataValue )
  {
    sum += 2 * ( x22 - x12 );
    weight += 2;
  }
  else if ( x12 == mInputNodataValue && x32 != mInputNodataValue && x22 != mInputNodataValue )
  {
    sum += 2 * ( x32 - x22 );
    weight += 2;
  }
 
  //third row
  if ( x33 != mInputNodataValue && x13 != mInputNodataValue ) //the normal case
  {
    sum += ( x33 - x13 );
    weight += 2;
  }
  else if ( x33 == mInputNodataValue && x13 != mInputNodataValue && x23 != mInputNodataValue )
  {
    sum += ( x23 - x13 );
    weight += 1;
  }
  else if ( x13 == mInputNodataValue && x33 != mInputNodataValue && x23 != mInputNodataValue )
  {
    sum += ( x33 - x23 );
    weight += 1;
  }
 
  if ( weight == 0 )
  {
    return mOutputNodataValue;
  }
 
  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 ) const
{
  //the basic formula would be simple, but we need to test for nodata values...
  //return (((*x11 - *x13) + 2 * (*x21 - *x23) + (*x31 - *x33)) / ( 8 * mCellSizeY));
 
  double sum = 0;
  int weight = 0;
 
  //first row
  if ( x11 != mInputNodataValue && x13 != mInputNodataValue ) //normal case
  {
    sum += ( x11 - x13 );
    weight += 2;
  }
  else if ( x11 == mInputNodataValue && x13 != mInputNodataValue && x12 != mInputNodataValue )
  {
    sum += ( x12 - x13 );
    weight += 1;
  }
  else if ( x31 == mInputNodataValue && x11 != mInputNodataValue && x12 != mInputNodataValue )
  {
    sum += ( x11 - x12 );
    weight += 1;
  }
 
  //second row
  if ( x21 != mInputNodataValue && x23 != mInputNodataValue )
  {
    sum += 2 * ( x21 - x23 );
    weight += 4;
  }
  else if ( x21 == mInputNodataValue && x23 != mInputNodataValue && x22 != mInputNodataValue )
  {
    sum += 2 * ( x22 - x23 );
    weight += 2;
  }
  else if ( x23 == mInputNodataValue && x21 != mInputNodataValue && x22 != mInputNodataValue )
  {
    sum += 2 * ( x21 - x22 );
    weight += 2;
  }
 
  //third row
  if ( x31 != mInputNodataValue && x33 != mInputNodataValue )
  {
    sum += ( x31 - x33 );
    weight += 2;
  }
  else if ( x31 == mInputNodataValue && x33 != mInputNodataValue && x32 != mInputNodataValue )
  {
    sum += ( x32 - x33 );
    weight += 1;
  }
  else if ( x33 == mInputNodataValue && x31 != mInputNodataValue && x32 != mInputNodataValue )
  {
    sum += ( x31 - x32 );
    weight += 1;
  }
 
  if ( weight == 0 )
  {
    return mOutputNodataValue;
  }
 
  return sum / ( weight * mCellSizeY ) * mZFactor;
}