108 auto outputBlock = std::make_unique<QgsRasterBlock>();
112 return outputBlock.release();
115 std::shared_ptr< QgsRasterBlock > inputBlock(
mInput->block( mBand,
extent, width, height, feedback ) );
117 if ( !inputBlock || inputBlock->isEmpty() )
120 return outputBlock.release();
123 std::shared_ptr< QgsRasterBlock > alphaBlock;
128 if ( !alphaBlock || alphaBlock->isEmpty() )
131 return outputBlock.release();
136 alphaBlock = inputBlock;
141 return outputBlock.release();
144 if ( width == 0 || height == 0 )
145 return outputBlock.release();
150 float cellXSize =
static_cast<float>(
extent.width() ) / width;
151 float cellYSize =
static_cast<float>(
extent.height() ) / height;
152 float zenithRad =
static_cast<float>( std::max( 0.0, 90 - mLightAngle ) * M_PI / 180.0 );
153 float azimuthRad =
static_cast<float>( -1 * mLightAzimuth * M_PI / 180.0 );
154 float cosZenithRad = std::cos( zenithRad );
155 float sinZenithRad = std::sin( zenithRad );
158 float cos_alt_mul_z = cosZenithRad *
static_cast<float>( mZFactor );
159 float cos_az_mul_cos_alt_mul_z = std::cos( azimuthRad ) * cos_alt_mul_z;
160 float sin_az_mul_cos_alt_mul_z = std::sin( azimuthRad ) * cos_alt_mul_z;
161 float cos_az_mul_cos_alt_mul_z_mul_254 = 254.0f * cos_az_mul_cos_alt_mul_z;
162 float sin_az_mul_cos_alt_mul_z_mul_254 = 254.0f * sin_az_mul_cos_alt_mul_z;
163 float square_z =
static_cast<float>( mZFactor * mZFactor );
164 float sin_altRadians_mul_254 = 254.0f * sinZenithRad;
167 float sin_altRadians_mul_127 = 127.0f * sinZenithRad;
169 float cos225_az_mul_cos_alt_mul_z_mul_127 = -32.87001872802012f * cos_alt_mul_z;
170 float cos_alt_mul_z_mul_127 = 127.0f * cos_alt_mul_z;
181 && inputBlock->dataTypeSize() <= 4 );
187 if ( source.isEmpty() )
195 std::chrono::time_point<std::chrono::system_clock> startTime( std::chrono::system_clock::now() );
203 std::size_t inputDataTypeSize = inputBlock->dataTypeSize();
204 std::size_t outputDataTypeSize = outputBlock->dataTypeSize();
207 switch ( inputBlock->dataType() )
210 typeName = u
"unsigned char"_s;
213 typeName = u
"unsigned int"_s;
216 typeName = u
"short"_s;
219 typeName = u
"unsigned int"_s;
225 typeName = u
"float"_s;
228 throw QgsException( u
"Unsupported data type for OpenCL processing."_s );
233 source.replace(
"__global float *scanLine"_L1, u
"__global %1 *scanLine"_s.arg( typeName ) );
237 std::size_t scanLineWidth( inputBlock->width() + 2 );
238 std::size_t inputSize( inputDataTypeSize * inputBlock->width() );
241 std::size_t bufferWidth( width + 2 );
242 std::size_t bufferSize( inputDataTypeSize * bufferWidth );
247 auto scanLine = std::make_unique<QgsRasterBlock>( inputBlock->dataType(), scanLineWidth, 1 );
254 std::vector<float> rasterParams;
255 rasterParams.push_back( inputBlock->noDataValue() );
256 rasterParams.push_back( outputBlock->noDataValue() );
257 rasterParams.push_back( mZFactor );
258 rasterParams.push_back( cellXSize );
259 rasterParams.push_back( cellYSize );
260 rasterParams.push_back(
static_cast<float>(
mOpacity ) );
263 rasterParams.push_back( cos_az_mul_cos_alt_mul_z_mul_254 );
264 rasterParams.push_back( sin_az_mul_cos_alt_mul_z_mul_254 );
265 rasterParams.push_back( square_z );
266 rasterParams.push_back( sin_altRadians_mul_254 );
269 rasterParams.push_back( sin_altRadians_mul_127 );
270 rasterParams.push_back( cos225_az_mul_cos_alt_mul_z_mul_127 );
271 rasterParams.push_back( cos_alt_mul_z_mul_127 );
274 rasterParams.push_back(
static_cast<float>( qBlue( defaultNodataColor ) ) );
275 rasterParams.push_back(
static_cast<float>( qGreen( defaultNodataColor ) ) );
276 rasterParams.push_back(
static_cast<float>( qRed( defaultNodataColor ) ) );
277 rasterParams.push_back(
static_cast<float>( qAlpha( defaultNodataColor ) ) / 255.0f );
280 rasterParams.push_back(
static_cast<float>( mMultiDirectional ) );
282 cl::Buffer rasterParamsBuffer( queue, rasterParams.begin(), rasterParams.end(),
true,
false,
nullptr );
283 cl::Buffer scanLine1Buffer( ctx, CL_MEM_READ_ONLY, bufferSize,
nullptr,
nullptr );
284 cl::Buffer scanLine2Buffer( ctx, CL_MEM_READ_ONLY, bufferSize,
nullptr,
nullptr );
285 cl::Buffer scanLine3Buffer( ctx, CL_MEM_READ_ONLY, bufferSize,
nullptr,
nullptr );
286 cl::Buffer *scanLineBuffer[3] = {&scanLine1Buffer, &scanLine2Buffer, &scanLine3Buffer};
288 cl::Buffer resultLineBuffer( ctx, CL_MEM_WRITE_ONLY, outputDataTypeSize * width,
nullptr,
nullptr );
290 static std::map<Qgis::DataType, cl::Program> programCache;
291 cl::Program program = programCache[inputBlock->dataType()];
292 if ( ! program.get() )
296 program = programCache[inputBlock->dataType()];
303 auto kernel = cl::KernelFunctor <
309 > ( program,
"processNineCellWindow" );
313 std::vector<int> rowIndex = {0, 1, 2};
315 for (
int i = 0; i < height; i++ )
324 feedback->
setProgress( 100.0 *
static_cast< double >( i ) / height );
330 scanLine->resetNoDataValue();
331 queue.enqueueWriteBuffer( scanLine1Buffer, CL_TRUE, 0, bufferSize, scanLine->bits( ) );
333 memcpy( scanLine->bits( 0, 1 ), inputBlock->bits( i, 0 ), inputSize );
334 queue.enqueueWriteBuffer( scanLine2Buffer, CL_TRUE, 0, bufferSize, scanLine->bits( ) );
336 memcpy( scanLine->bits( 0, 1 ), inputBlock->bits( i + 1, 0 ), inputSize );
337 queue.enqueueWriteBuffer( scanLine3Buffer, CL_TRUE, 0, bufferSize, scanLine->bits( ) );
343 if ( i == inputBlock->height() - 1 )
345 scanLine->resetNoDataValue();
346 queue.enqueueWriteBuffer( *scanLineBuffer[rowIndex[2]], CL_TRUE, 0, bufferSize, scanLine->bits( ) );
350 queue.enqueueWriteBuffer( *scanLineBuffer[rowIndex[2]], CL_TRUE, inputDataTypeSize * 1 , inputSize, inputBlock->bits( i + 1, 0 ) );
354 kernel( cl::EnqueueArgs(
358 *scanLineBuffer[rowIndex[0]],
359 *scanLineBuffer[rowIndex[1]],
360 *scanLineBuffer[rowIndex[2]],
365 queue.enqueueReadBuffer( resultLineBuffer, CL_TRUE, 0, outputDataTypeSize * outputBlock->width( ), outputBlock->bits( i, 0 ) );
366 std::rotate( rowIndex.begin(), rowIndex.begin() + 1, rowIndex.end() );
369 catch ( cl::Error &e )
384 double pixelValues[9] {0, 0, 0, 0, 0, 0, 0, 0, 0};
385 bool isNoData[9] {
false,
false,
false,
false,
false,
false,
false,
false,
false};
387 for (
int row = 0; row < height; row++ )
389 for (
int col = 0; col < width; col++ )
397 else if ( row == height - 1 )
405 pixelValues[ 0 ] = inputBlock->valueAndNoData( iUp, 0, isNoData[0] );
406 pixelValues[ 1 ] = pixelValues[0];
407 isNoData[1] = isNoData[0];
408 pixelValues[ 2 ] = pixelValues[0];
409 isNoData[2] = isNoData[0];
411 pixelValues[ 3 ] = inputBlock->valueAndNoData( row, 0, isNoData[3] );
412 pixelValues[ 4 ] = pixelValues[3];
413 isNoData[4] = isNoData[3];
414 pixelValues[ 5 ] = pixelValues[3];
415 isNoData[5] = isNoData[3];
417 pixelValues[ 6 ] = inputBlock->valueAndNoData( iDown, 0, isNoData[6] );
418 pixelValues[ 7 ] = pixelValues[6];
419 isNoData[7] = isNoData[6];
420 pixelValues[ 8 ] = pixelValues[6];
421 isNoData[8] = isNoData[6];
426 pixelValues[ 0 ] = pixelValues[1];
427 pixelValues[ 1 ] = pixelValues[2];
428 pixelValues[ 3 ] = pixelValues[4];
429 pixelValues[ 4 ] = pixelValues[5];
430 pixelValues[ 6 ] = pixelValues[7];
431 pixelValues[ 7 ] = pixelValues[8];
432 isNoData[ 0 ] = isNoData[1];
433 isNoData[ 1 ] = isNoData[2];
434 isNoData[ 3 ] = isNoData[4];
435 isNoData[ 4 ] = isNoData[5];
436 isNoData[ 6 ] = isNoData[7];
437 isNoData[ 7 ] = isNoData[8];
441 if ( col < width - 1 )
443 pixelValues[2] = inputBlock->valueAndNoData( iUp, col + 1, isNoData[2] );
444 pixelValues[5] = inputBlock->valueAndNoData( row, col + 1, isNoData[5] );
445 pixelValues[8] = inputBlock->valueAndNoData( iDown, col + 1, isNoData[8] );
450 outputBlock->setColor( row, col, defaultNodataColor );
455 const double x22 = pixelValues[4];
457 const double x11 = isNoData[0] ? x22 : pixelValues[0];
458 const double x21 = isNoData[3] ? x22 : pixelValues[3];
459 const double x31 = isNoData[6] ? x22 : pixelValues[6];
460 const double x12 = isNoData[1] ? x22 : pixelValues[1];
462 const double x32 = isNoData[7] ? x22 : pixelValues[7];
463 const double x13 = isNoData[2] ? x22 : pixelValues[2];
464 const double x23 = isNoData[5] ? x22 : pixelValues[5];
465 const double x33 = isNoData[8] ? x22 : pixelValues[8];
468 const double derX = ( ( x13 + x23 + x23 + x33 ) - ( x11 + x21 + x21 + x31 ) ) / ( 8 * cellXSize );
469 const double derY = ( ( x31 + x32 + x32 + x33 ) - ( x11 + x12 + x12 + x13 ) ) / ( 8 * -cellYSize );
474 if ( !mMultiDirectional )
477 grayValue = std::clamp( ( sin_altRadians_mul_254 -
478 ( derY * cos_az_mul_cos_alt_mul_z_mul_254 -
479 derX * sin_az_mul_cos_alt_mul_z_mul_254 ) ) /
480 std::sqrt( 1 + square_z * ( derX * derX + derY * derY ) ),
487 const float xx = derX * derX;
488 const float yy = derY * derY;
489 const float xx_plus_yy = xx + yy;
491 if ( xx_plus_yy == 0.0 )
493 grayValue = std::clamp(
static_cast<float>( 1.0 + sin_altRadians_mul_254 ), 0.0f, 255.0f );
498 float val225_mul_127 = sin_altRadians_mul_127 +
499 ( derX - derY ) * cos225_az_mul_cos_alt_mul_z_mul_127;
500 val225_mul_127 = ( val225_mul_127 <= 0.0 ) ? 0.0 : val225_mul_127;
501 float val270_mul_127 = sin_altRadians_mul_127 -
502 derX * cos_alt_mul_z_mul_127;
503 val270_mul_127 = ( val270_mul_127 <= 0.0 ) ? 0.0 : val270_mul_127;
504 float val315_mul_127 = sin_altRadians_mul_127 +
505 ( derX + derY ) * cos225_az_mul_cos_alt_mul_z_mul_127;
506 val315_mul_127 = ( val315_mul_127 <= 0.0 ) ? 0.0 : val315_mul_127;
507 float val360_mul_127 = sin_altRadians_mul_127 -
508 derY * cos_alt_mul_z_mul_127;
509 val360_mul_127 = ( val360_mul_127 <= 0.0 ) ? 0.0 : val360_mul_127;
512 const float weight_225 = 0.5 * xx_plus_yy - derX * derY;
513 const float weight_270 = xx;
514 const float weight_315 = xx_plus_yy - weight_225;
515 const float weight_360 = yy;
516 const float cang_mul_127 = (
517 ( weight_225 * val225_mul_127 +
518 weight_270 * val270_mul_127 +
519 weight_315 * val315_mul_127 +
520 weight_360 * val360_mul_127 ) / xx_plus_yy ) /
521 ( 1 + square_z * xx_plus_yy );
523 grayValue = std::clamp( 1.0f + cang_mul_127, 0.0f, 255.0f );
534 currentAlpha *= alphaBlock->value(
static_cast<qgssize>( row ) * width + col ) / 255.0;
539 outputBlock->setColor( row, col, qRgba( grayValue, grayValue, grayValue, 255 ) );
543 outputBlock->setColor( row, col, qRgba( currentAlpha * grayValue, currentAlpha * grayValue, currentAlpha * grayValue, currentAlpha * 255 ) );
552 if (
QgsSettings().value( u
"Map/logCanvasRefreshEvent"_s,
false ).toBool() )
555 .arg( useOpenCL ? u
"OpenCL"_s : u
"CPU"_s )
558 .arg( std::chrono::duration_cast<std::chrono::milliseconds>( std::chrono::system_clock::now() - startTime ).count() ),
564 return outputBlock.release();
616 QDomNodeList elements = element.elementsByTagName( u
"sld:RasterSymbolizer"_s );
617 if ( elements.size() == 0 )
621 QDomElement rasterSymbolizerElem = elements.at( 0 ).toElement();
628 QDomElement channelSelectionElem = doc.createElement( u
"sld:ChannelSelection"_s );
629 elements = rasterSymbolizerElem.elementsByTagName( u
"sld:Opacity"_s );
630 if ( elements.size() != 0 )
632 rasterSymbolizerElem.insertAfter( channelSelectionElem, elements.at( 0 ) );
636 elements = rasterSymbolizerElem.elementsByTagName( u
"sld:Geometry"_s );
637 if ( elements.size() != 0 )
639 rasterSymbolizerElem.insertAfter( channelSelectionElem, elements.at( 0 ) );
643 rasterSymbolizerElem.insertBefore( channelSelectionElem, rasterSymbolizerElem.firstChild() );
648 QDomElement channelElem = doc.createElement( u
"sld:GrayChannel"_s );
649 channelSelectionElem.appendChild( channelElem );
652 QDomElement sourceChannelNameElem = doc.createElement( u
"sld:SourceChannelName"_s );
653 sourceChannelNameElem.appendChild( doc.createTextNode( QString::number( mBand ) ) );
654 channelElem.appendChild( sourceChannelNameElem );
658 QDomElement shadedReliefElem = doc.createElement( u
"sld:ShadedRelief"_s );
659 rasterSymbolizerElem.appendChild( shadedReliefElem );
662 QDomElement brightnessOnlyElem = doc.createElement( u
"sld:BrightnessOnly"_s );
663 brightnessOnlyElem.appendChild( doc.createTextNode( u
"true"_s ) );
664 shadedReliefElem.appendChild( brightnessOnlyElem );
667 QDomElement reliefFactorElem = doc.createElement( u
"sld:ReliefFactor"_s );
668 reliefFactorElem.appendChild( doc.createTextNode( QString::number(
zFactor() ) ) );
669 shadedReliefElem.appendChild( reliefFactorElem );
672 QDomElement altitudeVendorOptionElem = doc.createElement( u
"sld:VendorOption"_s );
673 altitudeVendorOptionElem.setAttribute( u
"name"_s, u
"altitude"_s );
674 altitudeVendorOptionElem.appendChild( doc.createTextNode( QString::number(
altitude() ) ) );
675 shadedReliefElem.appendChild( altitudeVendorOptionElem );
678 QDomElement azimutVendorOptionElem = doc.createElement( u
"sld:VendorOption"_s );
679 azimutVendorOptionElem.setAttribute( u
"name"_s, u
"azimuth"_s );
680 azimutVendorOptionElem.appendChild( doc.createTextNode( QString::number(
azimuth() ) ) );
681 shadedReliefElem.appendChild( azimutVendorOptionElem );
684 QDomElement multidirectionalVendorOptionElem = doc.createElement( u
"sld:VendorOption"_s );
685 multidirectionalVendorOptionElem.setAttribute( u
"name"_s, u
"multidirectional"_s );
686 multidirectionalVendorOptionElem.appendChild( doc.createTextNode( QString::number(
multiDirectional() ) ) );
687 shadedReliefElem.appendChild( multidirectionalVendorOptionElem );