Skip to content
Permalink
Browse files

Merge pull request #7451 from elpaso/opencl-utils-2

[feature] OpenCL support
  • Loading branch information
elpaso committed Aug 8, 2018
2 parents 55473e7 + 0960b1f commit 0b502ff5b9b6925c694ae1709918e14233daf6cc
Showing with 2,752 additions and 331 deletions.
  1. +12 −0 CMakeLists.txt
  2. +2 −0 cmake_templates/qgsconfig.h.in
  3. +1 −0 images/images.qrc
  4. +14 −0 images/themes/default/mIconGPU.svg
  5. +1 −0 python/analysis/auto_generated/raster/qgsaspectfilter.sip.in
  6. +19 −3 python/analysis/auto_generated/raster/qgsninecellfilter.sip.in
  7. +2 −0 python/analysis/auto_generated/raster/qgsslopefilter.sip.in
  8. +12 −2 python/core/auto_generated/raster/qgsrasterblock.sip.in
  9. +3 −0 resources/CMakeLists.txt
  10. +40 −0 resources/opencl_programs/aspect.cl
  11. +53 −0 resources/opencl_programs/aspect_renderer.cl
  12. +71 −0 resources/opencl_programs/calcfirstder.cl
  13. +56 −0 resources/opencl_programs/hillshade.cl
  14. +110 −0 resources/opencl_programs/hillshade_renderer.cl
  15. +51 −0 resources/opencl_programs/ruggedness.cl
  16. +42 −0 resources/opencl_programs/slope.cl
  17. +45 −0 resources/opencl_programs/slope_renderer.cl
  18. +11 −0 src/analysis/CMakeLists.txt
  19. +1 −0 src/analysis/raster/qgsalignraster.cpp
  20. +10 −0 src/analysis/raster/qgsaspectfilter.h
  21. +36 −6 src/analysis/raster/qgshillshadefilter.cpp
  22. +23 −2 src/analysis/raster/qgshillshadefilter.h
  23. +346 −106 src/analysis/raster/qgsninecellfilter.cpp
  24. +55 −4 src/analysis/raster/qgsninecellfilter.h
  25. +0 −5 src/analysis/raster/qgsruggednessfilter.cpp
  26. +8 −0 src/analysis/raster/qgsruggednessfilter.h
  27. +4 −2 src/analysis/raster/qgsslopefilter.cpp
  28. +11 −0 src/analysis/raster/qgsslopefilter.h
  29. +35 −0 src/app/main.cpp
  30. +11 −0 src/app/qgisapp.cpp
  31. +49 −0 src/app/qgsoptions.cpp
  32. +25 −1 src/core/CMakeLists.txt
  33. +530 −0 src/core/qgsopenclutils.cpp
  34. +289 −0 src/core/qgsopenclutils.h
  35. +3 −4 src/core/raster/qgscolorrampshader.cpp
  36. +326 −138 src/core/raster/qgshillshaderenderer.cpp
  37. +12 −2 src/core/raster/qgsrasterblock.h
  38. +71 −32 src/ui/qgsoptionsbase.ui
  39. +112 −24 tests/src/analysis/testqgsninecellfilters.cpp
  40. +8 −0 tests/src/core/CMakeLists.txt
  41. +242 −0 tests/src/core/testqgsopenclutils.cpp
  42. BIN tests/testdata/analysis/hillshade.tif
  43. BIN tests/testdata/opencl/aspect.tif
  44. BIN tests/testdata/opencl/dem.tif
  45. BIN tests/testdata/opencl/slope.tif
@@ -28,6 +28,18 @@ MESSAGE(STATUS "QGIS version: ${COMPLETE_VERSION} ${RELEASE_NAME} (${QGIS_VERSIO


#############################################################
# Configure OpenCL if available

OPTION(USE_OPENCL "Use OpenCL" ON)
IF (USE_OPENCL)
FIND_PACKAGE(OpenCL)
IF(${OpenCL_FOUND})
SET(HAVE_OPENCL TRUE)
ELSE(${OpenCL_FOUND})
MESSAGE(STATUS "Couldn't find OpenCL: support DISABLED")
ENDIF(${OpenCL_FOUND})
ENDIF(USE_OPENCL)

# Configure CCache if available
IF(NOT MSVC)
option(USE_CCACHE "Use ccache" ON)
@@ -57,6 +57,8 @@

#cmakedefine HAVE_SERVER_PYTHON_PLUGINS

#cmakedefine HAVE_OPENCL

#cmakedefine ENABLE_MODELTEST

#cmakedefine HAVE_3D
@@ -699,6 +699,7 @@
<file>themes/default/mIndicatorEmbedded.svg</file>
<file>themes/default/mIconHistory.svg</file>
<file>themes/default/mIndicatorMemory.svg</file>
<file>themes/default/mIconGPU.svg</file>
</qresource>
<qresource prefix="/images/tips">
<file alias="symbol_levels.png">qgis_tips/symbol_levels.png</file>
@@ -0,0 +1,14 @@
<svg xmlns="http://www.w3.org/2000/svg" viewBox="0 0 22 22">
<defs id="defs3051">
<style type="text/css" id="current-color-scheme">
.ColorScheme-Text {
color:#4d4d4d;
}
</style>
</defs>
<path
style="fill:currentColor;fill-opacity:1;stroke:none"
d="M 6 3 L 6 5 L 5 5 L 5 6 L 3 6 L 3 7 L 5 7 L 5 9 L 3 9 L 3 10 L 5 10 L 5 12 L 3 12 L 3 13 L 5 13 L 5 15 L 3 15 L 3 16 L 5 16 L 5 17 L 6 17 L 6 19 L 7 19 L 7 17 L 9 17 L 9 19 L 10 19 L 10 17 L 12 17 L 12 19 L 13 19 L 13 17 L 15 17 L 15 19 L 16 19 L 16 17 L 17 17 L 17 16 L 19 16 L 19 15 L 17 15 L 17 13 L 19 13 L 19 12 L 17 12 L 17 10 L 19 10 L 19 9 L 17 9 L 17 7 L 19 7 L 19 6 L 17 6 L 17 5 L 16 5 L 16 3 L 15 3 L 15 5 L 13 5 L 13 3 L 12 3 L 12 5 L 10 5 L 10 3 L 9 3 L 9 5 L 7 5 L 7 3 L 6 3 z "
class="ColorScheme-Text"
/>
</svg>
@@ -29,6 +29,7 @@ Calculates output value from nine input values. The input values and the output
nodata value if not present or outside of the border. Must be implemented by subclasses*
%End


};

/************************************************************************
@@ -10,6 +10,7 @@




class QgsNineCellFilter
{
%Docstring
@@ -34,7 +35,7 @@ Starts the calculation, reads from mInputFile and stores the result in mOutputFi

:param feedback: feedback object that receives update and that is checked for cancelation.

:return: 0 in case of success*
:return: 0 in case of success
%End

double cellSizeX() const;
@@ -54,8 +55,23 @@ Starts the calculation, reads from mInputFile and stores the result in mOutputFi
float *x12, float *x22, float *x32,
float *x13, float *x23, float *x33 ) = 0;
%Docstring
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*
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.

First index of the input cell is the row, second index is the column

:param x11: surrounding cell top left
:param x21: surrounding cell central left
:param x31: surrounding cell bottom left
:param x12: surrounding cell top central
:param x22: the central cell for which the value will be calculated
:param x32: surrounding cell bottom central
:param x13: surrounding cell top right
:param x23: surrounding cell central right
:param x33: surrounding cell bottom right

:return: the calculated cell value for the central cell x22
%End

protected:
@@ -28,6 +28,8 @@ Calculates slope values in a window of 3x3 cells based on first order derivative
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*
%End


};

/************************************************************************
@@ -186,7 +186,7 @@ Read a single value
:return: color *
%End

bool isNoData( int row, int column );
bool isNoData( int row, int column ) const;
%Docstring
Check if value at position is no data

@@ -196,7 +196,17 @@ Check if value at position is no data
:return: true if value is no data *
%End

bool isNoData( qgssize index );
bool isNoData( qgssize row, qgssize column ) const;
%Docstring
Check if value at position is no data

:param row: row index
:param column: column index

:return: true if value is no data *
%End

bool isNoData( qgssize index ) const;
%Docstring
Check if value at position is no data

@@ -11,6 +11,9 @@ INSTALL(DIRECTORY themes DESTINATION ${QGIS_DATA_DIR}/resources)
INSTALL(DIRECTORY data DESTINATION ${QGIS_DATA_DIR}/resources)
INSTALL(DIRECTORY metadata-ISO DESTINATION ${QGIS_DATA_DIR}/resources)
INSTALL(DIRECTORY palettes DESTINATION ${QGIS_DATA_DIR}/resources)
IF (HAVE_OPENCL)
INSTALL(DIRECTORY opencl_programs DESTINATION ${QGIS_DATA_DIR}/resources)
ENDIF (HAVE_OPENCL)

IF (WITH_SERVER)
INSTALL(DIRECTORY server DESTINATION ${QGIS_DATA_DIR}/resources)
@@ -0,0 +1,40 @@
#include "calcfirstder.cl"

__kernel void processNineCellWindow( __global float *scanLine1,
__global float *scanLine2,
__global float *scanLine3,
__global float *resultLine,
__global float *rasterParams // mInputNodataValue, mOutputNodataValue, mZFactor, mCellSizeX, mCellSizeY
)
{

// Get the index of the current element
const int i = get_global_id(0);

if ( scanLine2[i+1] == rasterParams[0] )
{
resultLine[i] = rasterParams[1];
}
else
{
float derX = calcFirstDer( scanLine1[i], scanLine1[i+1], scanLine1[i+2],
scanLine2[i], scanLine2[i+1], scanLine2[i+2],
scanLine3[i], scanLine3[i+1], scanLine3[i+2],
rasterParams[0], rasterParams[1], rasterParams[2], rasterParams[3] );

float derY = calcFirstDer( scanLine3[i], scanLine2[i], scanLine1[i],
scanLine3[i+1], scanLine2[i+1], scanLine1[i+1],
scanLine3[i+2], scanLine2[i+2], scanLine1[i+2],
rasterParams[0], rasterParams[1], rasterParams[2], rasterParams[4]);

if ( derX == rasterParams[1] || derY == rasterParams[1] ||
( derX == 0.0f && derY == 0.0f) )
{
resultLine[i] = rasterParams[1];
}
else
{
resultLine[i] = 180.0f + atan2pi( derX, derY ) * 180.0f;
}
}
}
@@ -0,0 +1,53 @@
#include "calcfirstder.cl"

// Aspect renderer for QGIS

__kernel void processNineCellWindow( __global float *scanLine1,
__global float *scanLine2,
__global float *scanLine3,
__global uchar4 *resultLine,
__global float *rasterParams
) {

// Get the index of the current element
const int i = get_global_id(0);

// Do the operation
//return (( (x31 - x11) + 2 * (x32 - x12) + (x33 - x13) ) / (8 * mCellSizeX))
float derX = calcFirstDer( scanLine1[i], scanLine2[i], scanLine3[i],
scanLine1[i+1], scanLine2[i+1], scanLine3[i+1],
scanLine1[i+2], scanLine2[i+2], scanLine3[i+2],
rasterParams[0], rasterParams[1], rasterParams[2], rasterParams[3]
);
//return (((x11 - x13) + 2 * (x21 - x23) + (x31 - x33)) / ( 8 * mCellSizeY));
float derY = calcFirstDer( scanLine1[i+2], scanLine1[i+1], scanLine1[i],
scanLine2[i+2], scanLine2[i+1], scanLine2[i],
scanLine3[i+2], scanLine3[i+1], scanLine3[i],
rasterParams[0], rasterParams[1], rasterParams[2], rasterParams[4]
);


float res;
if ( derX == rasterParams[1] || derY == rasterParams[1] ||
( derX == 0.0f && derY == 0.0f) )
{
res = rasterParams[1];
}
else
{
// 180.0 / M_PI = 57.29577951308232
float aspect = atan2( derX, derY ) * 57.29577951308232;
if ( aspect < 0 )
res = 90.0f - aspect;
else if (aspect > 90.0f)
// 360 + 90 = 450
res = 450.0f - aspect;
else
res = 90.0 - aspect;
}

res = res / 360 * 255;

resultLine[i] = (uchar4)(res, res, res, 255);
}

@@ -0,0 +1,71 @@
// Calculate the first derivative from a 3x3 cell matrix
float calcFirstDer( float x11, float x21, float x31, float x12, float x22, float x32, float x13, float x23, float x33,
float inputNodataValue, float outputNodataValue, float zFactor, float mCellSize )
{
//the basic formula would be simple, but we need to test for nodata values...
//X: return (( (x31 - x11) + 2 * (x32 - x12) + (x33 - x13) ) / (8 * cellSizeX));
//Y: return (((x11 - x13) + 2 * (x21 - x23) + (x31 - x33)) / ( 8 * cellSizeY));

int weight = 0;
float sum = 0;


//first row
if ( x31 != inputNodataValue && x11 != inputNodataValue ) //the normal case
{
sum += ( x31 - x11 );
weight += 2;
}
else if ( x31 == inputNodataValue && x11 != inputNodataValue && x21 != inputNodataValue ) //probably 3x3 window is at the border
{
sum += ( x21 - x11 );
weight += 1;
}
else if ( x11 == inputNodataValue && x31 != inputNodataValue && x21 != inputNodataValue ) //probably 3x3 window is at the border
{
sum += ( x31 - x21 );
weight += 1;
}

//second row
if ( x32 != inputNodataValue && x12 != inputNodataValue ) //the normal case
{
sum += 2.0f * ( x32 - x12 );
weight += 4;
}
else if ( x32 == inputNodataValue && x12 != inputNodataValue && x22 != inputNodataValue )
{
sum += 2.0f * ( x22 - x12 );
weight += 2;
}
else if ( x12 == inputNodataValue && x32 != inputNodataValue && x22 != inputNodataValue )
{
sum += 2.0f * ( x32 - x22 );
weight += 2;
}

//third row
if ( x33 != inputNodataValue && x13 != inputNodataValue ) //the normal case
{
sum += ( x33 - x13 );
weight += 2;
}
else if ( x33 == inputNodataValue && x13 != inputNodataValue && x23 != inputNodataValue )
{
sum += ( x23 - x13 );
weight += 1;
}
else if ( x13 == inputNodataValue && x33 != inputNodataValue && x23 != inputNodataValue )
{
sum += ( x33 - x23 );
weight += 1;
}

if ( weight == 0 )
{
return outputNodataValue;
}

return sum / ( weight * mCellSize ) * zFactor;
}

@@ -0,0 +1,56 @@
#include "calcfirstder.cl"

#ifndef M_PI
#define M_PI 3.14159265358979323846
#endif

__kernel void processNineCellWindow( __global float *scanLine1,
__global float *scanLine2,
__global float *scanLine3,
__global float *resultLine,
__global float *rasterParams // mInputNodataValue, mOutputNodataValue, mZFactor, mCellSizeX, mCellSizeY, zenith_rad, azimuth_rad
)
{

// Get the index of the current element
const int i = get_global_id(0);

if ( scanLine2[i+1] == rasterParams[0] )
{
resultLine[i] = rasterParams[1];
}
else
{

float derX = calcFirstDer( scanLine1[i], scanLine1[i+1], scanLine1[i+2],
scanLine2[i], scanLine2[i+1], scanLine2[i+2],
scanLine3[i], scanLine3[i+1], scanLine3[i+2],
rasterParams[0], rasterParams[1], rasterParams[2], rasterParams[3] );

float derY = calcFirstDer( scanLine3[i], scanLine2[i], scanLine1[i],
scanLine3[i+1], scanLine2[i+1], scanLine1[i+1],
scanLine3[i+2], scanLine2[i+2], scanLine1[i+2],
rasterParams[0], rasterParams[1], rasterParams[2], rasterParams[4]);

if ( derX == rasterParams[1] || derY == rasterParams[1] )
{
resultLine[i] = rasterParams[1];
}
else
{

float slope_rad = sqrt( derX * derX + derY * derY );
slope_rad = atan( slope_rad );
float aspect_rad;
if ( derX == 0.0f && derY == 0.0f)
aspect_rad = rasterParams[7] / 2.0f;
else
aspect_rad = M_PI + atan2( derX, derY );

resultLine[i] = max(0.0f, 255.0f * ( ( rasterParams[5] * cos( slope_rad ) ) +
( rasterParams[6] * sin( slope_rad ) *
cos( rasterParams[7] - aspect_rad ) ) ) );

}
}
}

0 comments on commit 0b502ff

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