| 240 | */ |
| 241 | |
| 242 | CPLErr GDALGridInverseDistanceToAPowerNearestNeighbor( |
| 243 | const void *poOptionsIn, GUInt32 nPoints, const double *padfX, |
| 244 | const double *padfY, const double *padfZ, double dfXPoint, double dfYPoint, |
| 245 | double *pdfValue, void *hExtraParamsIn) |
| 246 | { |
| 247 | CPL_IGNORE_RET_VAL(nPoints); |
| 248 | |
| 249 | const GDALGridInverseDistanceToAPowerNearestNeighborOptions *const |
| 250 | poOptions = static_cast< |
| 251 | const GDALGridInverseDistanceToAPowerNearestNeighborOptions *>( |
| 252 | poOptionsIn); |
| 253 | const double dfRadius = poOptions->dfRadius; |
| 254 | const double dfSmoothing = poOptions->dfSmoothing; |
| 255 | const double dfSmoothing2 = dfSmoothing * dfSmoothing; |
| 256 | |
| 257 | const GUInt32 nMaxPoints = poOptions->nMaxPoints; |
| 258 | |
| 259 | GDALGridExtraParameters *psExtraParams = |
| 260 | static_cast<GDALGridExtraParameters *>(hExtraParamsIn); |
| 261 | const CPLQuadTree *phQuadTree = psExtraParams->hQuadTree; |
| 262 | CPLAssert(phQuadTree); |
| 263 | |
| 264 | const double dfRPower2 = psExtraParams->dfRadiusPower2PreComp; |
| 265 | const double dfPowerDiv2 = psExtraParams->dfPowerDiv2PreComp; |
| 266 | |
| 267 | std::multimap<double, double> oMapDistanceToZValues; |
| 268 | |
| 269 | const double dfSearchRadius = dfRadius; |
| 270 | CPLRectObj sAoi; |
| 271 | sAoi.minx = dfXPoint - dfSearchRadius; |
| 272 | sAoi.miny = dfYPoint - dfSearchRadius; |
| 273 | sAoi.maxx = dfXPoint + dfSearchRadius; |
| 274 | sAoi.maxy = dfYPoint + dfSearchRadius; |
| 275 | int nFeatureCount = 0; |
| 276 | GDALGridPoint **papsPoints = reinterpret_cast<GDALGridPoint **>( |
| 277 | CPLQuadTreeSearch(phQuadTree, &sAoi, &nFeatureCount)); |
| 278 | if (nFeatureCount != 0) |
| 279 | { |
| 280 | for (int k = 0; k < nFeatureCount; k++) |
| 281 | { |
| 282 | const int i = papsPoints[k]->i; |
| 283 | const double dfRX = padfX[i] - dfXPoint; |
| 284 | const double dfRY = padfY[i] - dfYPoint; |
| 285 | |
| 286 | const double dfR2 = dfRX * dfRX + dfRY * dfRY; |
| 287 | // real distance + smoothing |
| 288 | const double dfRsmoothed2 = dfR2 + dfSmoothing2; |
| 289 | if (dfRsmoothed2 < 0.0000000000001) |
| 290 | { |
| 291 | *pdfValue = padfZ[i]; |
| 292 | CPLFree(papsPoints); |
| 293 | return CE_None; |
| 294 | } |
| 295 | // is point within real distance? |
| 296 | if (dfR2 <= dfRPower2) |
| 297 | { |
| 298 | oMapDistanceToZValues.insert( |
| 299 | std::make_pair(dfRsmoothed2, padfZ[i])); |
nothing calls this directly
no test coverage detected