3

I have a large point vector layer (a few million objects) with measured elevation values and have converted it to raster, using the raster>conversion>(vector to raster) option in QGIS. It seems to work fine, but I haven't been able to find info on how it exactly the value for a pixel is computed.

Does QGIS take the average of all points in the cell, the median, the first object it finds in the attribute table, etc.?

Is there perhaps a plugin with options for this?

  • 1
    The joy of open-source is that you don't have to guess, just dig down deep enough in the source tree to locate the code. – Vince May 12 at 19:20
  • I suspect you'll find the scanline fill algorithm alive and well. – Vince May 12 at 21:37
  • @Vince the gvBurnScanLine function is definitely still there but the function used is GDALdllImageLine or GDALdllImageLineAllTouched that is actually used in gv_rasterize_one_shape. – Michael Stimson May 13 at 1:14
6

In gdalrasterize.cpp, function CPLErr GDALRasterizeGeometries, line 715 to 742 in my local snapshot (line 885 to 914 of the online source at the time of writing):

    eErr =
        poDS->RasterIO(GF_Read,
                       0, iY, poDS->GetRasterXSize(), nThisYChunkSize,
                       pabyChunkBuf,
                       poDS->GetRasterXSize(), nThisYChunkSize,
                       eType, nBandCount, panBandList,
                       0, 0, 0, NULL);
    if( eErr != CE_None )
        break;

    for( int iShape = 0; iShape < nGeomCount; iShape++ )
    {
        gv_rasterize_one_shape( pabyChunkBuf, iY,
                                poDS->GetRasterXSize(), nThisYChunkSize,
                                nBandCount, eType, bAllTouched,
                                reinterpret_cast<OGRGeometry *>(
                                                    pahGeometries[iShape]),
                                padfGeomBurnValue + iShape*nBandCount,
                                eBurnValueSource, eMergeAlg,
                                pfnTransformer, pTransformArg );
    }

    eErr =
        poDS->RasterIO( GF_Write, 0, iY,
                        poDS->GetRasterXSize(), nThisYChunkSize,
                        pabyChunkBuf,
                        poDS->GetRasterXSize(), nThisYChunkSize,
                        eType, nBandCount, panBandList, 0, 0, 0, NULL);

For those that don't speak C the key is for( int iShape = 0; iShape < nGeomCount; iShape++ ) which means that values are overwritten without checking if the cell already has a value... so it's the last geometry in the cell that sets the final value.

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.