Here's what I've settled on.  Let me know your thoughts.  The algorithm now
correctly calculates angle to lowest nearby point.  If it's a local minima,
it instead calculates angle to the nearest high point, which will be a
negative number, so if the user wants to nullify these or use them in
additional processing steps they can.  How does adding code to gdal work?

Best,
Steve

float GDALSlopeHydroAlg (float* afWin, float fDstNoDataValue, void* pData)

{
    // Hydrologic Slope is the max
    //  local slope btw center cell and adjacent cells

    const double radiansToDegrees = 180.0 / M_PI;
    GDALSlopeAlgData* psData = (GDALSlopeAlgData*)pData;
                double key;

    float pafLocalMin = afWin[0];
     float pafLocalMax = afWin[0];

    for ( int k = 1; k < 9; k++)
    {
        if (afWin[k] < pafLocalMin)
        {
            pafLocalMin=afWin[k];
        }

        if (afWin[k] > pafLocalMax)
        {
            pafLocalMax=afWin[k];
        }

    }

    key = afWin[4] - pafLocalMin;

    if (key < 0)
     {
        key = afWin[4] - pafLocalMax;
     }

    if (psData->slopeFormat == 1)
        return (float) (atan(sqrt(key) / (2*psData->scale)) *
radiansToDegrees);
    else
        return (float) (100*(sqrt(key) / (2*psData->scale)));
}
_______________________________________________
gdal-dev mailing list
gdal-dev@lists.osgeo.org
http://lists.osgeo.org/mailman/listinfo/gdal-dev

Reply via email to