Hillshaderenderer: Add 'multidirectional oblique-weighted' rendering …

…option (#3160)

python/core/raster/qgshillshaderenderer.sip

@@ -46,6 +46,12 @@ class QgsHillshadeRenderer : QgsRasterRenderer
      */
     double zFactor()  const;
 
+    /**
+     * @brief Is Multi Directional
+     * @return Is Multi Directional
+     */
+    bool multiDirectional() const;
+
 
     /**
      * @brief Set the azimith of the light source.
@@ -64,4 +70,10 @@ class QgsHillshadeRenderer : QgsRasterRenderer
      * @param zfactor The z factor.
      */
     void setZFactor( double zfactor );
+
+    /**
+     * @brief Set Is Multi Directional
+     * @param isMultiDirectional Is multi directional
+     */
+    void setMultiDirectional( bool isMultiDirectional );
 };

python/gui/raster/qgshillshaderendererwidget.sip

@@ -51,6 +51,12 @@ class QgsHillshadeRendererWidget: QgsRasterRendererWidget
      */
     double zFactor()  const;
 
+    /**
+     * @brief multiDirectional
+     * @return multiDirectional
+     */
+    bool multiDirectional() const;
+
   public slots:
     /**
      * @brief Set the altitude of the light source
@@ -69,4 +75,10 @@ class QgsHillshadeRendererWidget: QgsRasterRendererWidget
      * @param zfactor The z factor.
      */
     void setZFactor( double zfactor );
+
+    /**
+     * @brief set MultiDirectional
+     * @param isMultiDirectional
+     */
+    void setMultiDirectional( bool isMultiDirectional);
 };

src/core/raster/qgshillshaderenderer.cpp

@@ -30,6 +30,7 @@ QgsHillshadeRenderer::QgsHillshadeRenderer( QgsRasterInterface *input, int band,
     , mZFactor( 1 )
     , mLightAngle( lightAngle )
     , mLightAzimuth( lightAzimuth )
+    , mMultiDirectional( false )
 {
 
 }
@@ -38,6 +39,7 @@ QgsHillshadeRenderer *QgsHillshadeRenderer::clone() const
 {
   QgsHillshadeRenderer* r = new QgsHillshadeRenderer( nullptr, mBand, mLightAzimuth, mLightAngle );
   r->setZFactor( mZFactor );
+  r->setMultiDirectional( mMultiDirectional );
   return r;
 }
 
@@ -52,8 +54,10 @@ QgsRasterRenderer *QgsHillshadeRenderer::create( const QDomElement &elem, QgsRas
   double azimuth = elem.attribute( "azimuth", "315" ).toDouble();
   double angle = elem.attribute( "angle", "45" ).toDouble();
   double zFactor = elem.attribute( "zfactor", "1" ).toDouble();
+  int multiDirection = elem.attribute( "multidirection", "0" ).toInt();
   QgsHillshadeRenderer* r = new QgsHillshadeRenderer( input, band, azimuth , angle );
   r->setZFactor( zFactor );
+  r->setMultiDirectional( multiDirection > 0 );
   return r;
 }
 
@@ -71,6 +75,7 @@ void QgsHillshadeRenderer::writeXML( QDomDocument &doc, QDomElement &parentElem
   rasterRendererElem.setAttribute( "azimuth", QString::number( mLightAzimuth ) );
   rasterRendererElem.setAttribute( "angle", QString::number( mLightAngle ) );
   rasterRendererElem.setAttribute( "zfactor", QString::number( mZFactor ) );
+  rasterRendererElem.setAttribute( "multidirection", QString::number( mMultiDirectional ) );
   parentElem.appendChild( rasterRendererElem );
 }
 
@@ -106,6 +111,12 @@ QgsRasterBlock *QgsHillshadeRenderer::block( int bandNo, const QgsRectangle &ext
   double cosZenithRad = cos( zenithRad );
   double sinZenithRad = sin( zenithRad );
 
+  // Multi direction hillshade: http://pubs.usgs.gov/of/1992/of92-422/of92-422.pdf
+  double angle0Rad = (-1 * mLightAzimuth - 45 - 45 * 0.5) * M_PI / 180.0;
+  double angle1Rad = (-1 * mLightAzimuth - 45 * 0.5) * M_PI / 180.0;
+  double angle2Rad = (-1 * mLightAzimuth + 45 * 0.5) * M_PI / 180.0;
+  double angle3Rad = (-1 * mLightAzimuth + 45 + 45 * 0.5) * M_PI / 180.0;
+
   QRgb myDefaultColor = NODATA_COLOR;
 
   for ( qgssize i = 0; i < ( qgssize )height; i++ )
@@ -181,14 +192,39 @@ QgsRasterBlock *QgsHillshadeRenderer::block( int bandNo, const QgsRectangle &ext
       double derX = calcFirstDerX( x11, x21, x31, x12, x22, x32, x13, x23, x33, cellXSize );
       double derY = calcFirstDerY( x11, x21, x31, x12, x22, x32, x13, x23, x33, cellYSize );
 
-      double slope_rad = atan( mZFactor * sqrt( derX * derX + derY * derY ) );
+      double slopeRad = atan( mZFactor * sqrt( derX * derX + derY * derY ) );
       double aspectRad = atan2( derX, -derY );
 
-      double colorvalue = qBound( 0.0, 255.0 * (( cosZenithRad * cos( slope_rad ) ) +
-                                  ( sinZenithRad * sin( slope_rad ) *
-                                    cos( azimuthRad - aspectRad ) ) ), 255.0 );
 
-      outputBlock->setColor( i, j, qRgb( colorvalue, colorvalue, colorvalue ) );
+      double grayValue;
+      if( !mMultiDirectional )
+      {
+          // Standard single direction hillshade
+          grayValue = qBound( 0.0, 255.0 * ( cosZenithRad * cos( slopeRad )
+                                             + sinZenithRad * sin( slopeRad )
+                                             * cos( azimuthRad - aspectRad )) , 255.0 );
+      }
+      else
+      {
+          // Weighted multi direction as in http://pubs.usgs.gov/of/1992/of92-422/of92-422.pdf
+          double weight0 = sin( aspectRad - angle0Rad );
+          double weight1 = sin( aspectRad - angle1Rad );
+          double weight2 = sin( aspectRad - angle2Rad );
+          double weight3 = sin( aspectRad - angle3Rad );
+          weight0 = weight0 * weight0;
+          weight1 = weight1 * weight1;
+          weight2 = weight2 * weight2;
+          weight3 = weight3 * weight3;
+
+          double cosSlope = cosZenithRad * cos( slopeRad );
+          double sinSlope = sinZenithRad * sin( slopeRad );
+          double color0 = cosSlope + sinSlope * cos( angle0Rad - aspectRad ) ;
+          double color1 = cosSlope + sinSlope * cos( angle1Rad - aspectRad ) ;
+          double color2 = cosSlope + sinSlope * cos( angle2Rad - aspectRad ) ;
+          double color3 = cosSlope + sinSlope * cos( angle3Rad - aspectRad ) ;
+          grayValue = qBound(0.0, 255 * ( weight0 * color0 + weight1 * color1 + weight2 * color2 + weight3 * color3 ) * 0.5, 255.0);
+      }
+      outputBlock->setColor( i, j, qRgb( grayValue, grayValue, grayValue ) );
     }
   }
   return outputBlock;

src/core/raster/qgshillshaderenderer.h

@@ -84,6 +84,12 @@ class CORE_EXPORT QgsHillshadeRenderer : public QgsRasterRenderer
      */
     double zFactor()  const { return mZFactor; }
 
+    /**
+     * @brief Is Multi Directional
+     * @return Is Multi Directional
+     */
+    bool multiDirectional() const { return mMultiDirectional; }
+
     /**
      * @brief Set the azimith of the light source.
      * @param azimuth The azimuth of the light source.
@@ -102,11 +108,18 @@ class CORE_EXPORT QgsHillshadeRenderer : public QgsRasterRenderer
      */
     void setZFactor( double zfactor ) { mZFactor = zfactor; }
 
+    /**
+     * @brief Set Is Multi Directional
+     * @param isMultiDirectional Is multi directional
+     */
+    void setMultiDirectional( bool isMultiDirectional ){ mMultiDirectional = isMultiDirectional; }
+
   private:
     int mBand;
     double mZFactor;
     double mLightAngle;
     double mLightAzimuth;
+    bool   mMultiDirectional;
 
     /** Calculates the first order derivative in x-direction according to Horn (1981)*/
     double calcFirstDerX( double x11, double x21, double x31, double x12, double x22, double x32, double x13, double x23, double x33 , double cellsize );

src/gui/raster/qgshillshaderendererwidget.cpp

@@ -40,6 +40,8 @@ QgsHillshadeRendererWidget::QgsHillshadeRendererWidget( QgsRasterLayer *layer, c
   mZFactor->setValue( 1 );
   mZFactor->setClearValue( 1 );
 
+  mMultiDirection->setChecked( false );
+
   if ( mRasterLayer )
   {
     QgsRasterDataProvider* provider = mRasterLayer->dataProvider();
@@ -63,6 +65,7 @@ QgsHillshadeRendererWidget::QgsHillshadeRendererWidget( QgsRasterLayer *layer, c
   connect( mLightAzimuth, SIGNAL( valueChanged( double ) ), this, SLOT( on_mLightAzimuth_updated( double ) ) );
   connect( mLightAzimuthDial, SIGNAL( valueChanged( int ) ), this, SLOT( on_mLightAzimuthDail_updated( int ) ) );
   connect( mZFactor, SIGNAL( valueChanged( double ) ), this, SIGNAL( widgetChanged() ) );
+  connect( mMultiDirection, SIGNAL( stateChanged( int )), this, SIGNAL( widgetChanged() ) );
 
   QgsBilinearRasterResampler* zoomedInResampler = new QgsBilinearRasterResampler();
   layer->resampleFilter()->setZoomedInResampler( zoomedInResampler );
@@ -91,6 +94,7 @@ QgsRasterRenderer *QgsHillshadeRendererWidget::renderer()
   QgsHillshadeRenderer* renderer = new QgsHillshadeRenderer( provider, band, mLightAzimuth->value(), mLightAngle->value() );
   double value = mZFactor->value();
   renderer->setZFactor( value );
+  renderer->setMultiDirectional( mMultiDirection->checkState() );
   return renderer;
 }
 
@@ -103,6 +107,7 @@ void QgsHillshadeRendererWidget::setFromRenderer( const QgsRasterRenderer *rende
     mLightAngle->setValue( r->altitude() );
     mLightAzimuth->setValue( r->azimuth() );
     mZFactor->setValue( r->zFactor() );
+    mMultiDirection->setChecked( r->multiDirectional() );
   }
 }
 
@@ -121,6 +126,11 @@ void QgsHillshadeRendererWidget::setZFactor( double zfactor )
   mZFactor->setValue( zfactor );
 }
 
+void QgsHillshadeRendererWidget::setMultiDirectional(bool isMultiDirectional )
+{
+    mMultiDirection->setChecked( isMultiDirectional );
+}
+
 void QgsHillshadeRendererWidget::on_mLightAzimuth_updated( double value )
 {
   int newvalue = ( int )value - 180;
@@ -153,3 +163,8 @@ double QgsHillshadeRendererWidget::zFactor() const
 {
   return mZFactor->value();
 }
+
+bool QgsHillshadeRendererWidget::multiDirectional() const
+{
+    return mMultiDirection->isChecked();
+}

src/gui/raster/qgshillshaderendererwidget.h

@@ -77,6 +77,12 @@ class GUI_EXPORT QgsHillshadeRendererWidget: public QgsRasterRendererWidget, pri
      */
     double zFactor()  const;
 
+    /**
+     * @brief multiDirection
+     * @return multiDirection
+     */
+    bool multiDirectional() const;
+
   public slots:
     /**
      * @brief Set the altitude of the light source
@@ -96,6 +102,12 @@ class GUI_EXPORT QgsHillshadeRendererWidget: public QgsRasterRendererWidget, pri
      */
     void setZFactor( double zfactor );
 
+    /**
+     * @brief set MultiDirection
+     * @param isMultiDirection
+     */
+    void setMultiDirectional( bool isMultiDirectional);
+
   private slots:
     void on_mLightAzimuth_updated( double value );
     void on_mLightAzimuthDail_updated( int value );

src/ui/raster/qghillshaderendererwidget.ui

@@ -99,7 +99,7 @@
      </property>
     </widget>
    </item>
-   <item row="4" column="0" colspan="3">
+   <item row="5" column="0" colspan="3">
     <spacer name="verticalSpacer">
      <property name="orientation">
       <enum>Qt::Vertical</enum>
@@ -135,6 +135,20 @@
      </property>
     </widget>
    </item>
+   <item row="4" column="2">
+    <widget class="QCheckBox" name="mMultiDirection">
+     <property name="text">
+      <string/>
+     </property>
+    </widget>
+   </item>
+   <item row="4" column="0">
+    <widget class="QLabel" name="label_5">
+     <property name="text">
+      <string>Multidirectional</string>
+     </property>
+    </widget>
+   </item>
   </layout>
  </widget>
  <customwidgets>