| 523 | } |
| 524 | |
| 525 | inline void GoldAddFeatures(RawBufferType &featCoords, const long2 &featCoordsStrides, RawBufferType &featMetadata, |
| 526 | const long2 &featMetadataStrides, RawBufferType &featDescriptors, |
| 527 | const long2 &featDescriptorsStrides, int maxCapacity, RawBufferType &numFeatures, |
| 528 | const long1 &numFeaturesStrides, const RawPyramidType &srcGaussianPyramid, |
| 529 | const RawPyramidType &srcDoGPyramid, const long3 &currStrides, const long3 &currShape, |
| 530 | int octave, int firstOctave, int numOctaveLayers, float contrastThreshold, |
| 531 | float edgeThreshold, float initSigma, int l, int currBatch, int r, int c) |
| 532 | { |
| 533 | constexpr float kImageScale = 1.f / cuda::TypeTraits<VT>::max; // source images data type scale |
| 534 | constexpr float kDScale1 = kImageScale * .5f; // first derivative scale |
| 535 | constexpr float kDScale2 = kImageScale; // second derivative scale |
| 536 | constexpr float kDScaleC = kImageScale * .25f; // cross derivative scale |
| 537 | |
| 538 | float cv; // central value |
| 539 | cuda::math::Vector<float, 3> dD, sol; // derivative distances and solver solution |
| 540 | cuda::math::Matrix<float, 3, 3> H; // Hessian matrix |
| 541 | |
| 542 | auto dogVal = [&srcDoGPyramid, &currStrides, &octave, &currBatch](int layer, int row, int col) |
| 543 | { |
| 544 | return util::ValueAt<WT>(srcDoGPyramid[octave][layer], currStrides, long3{currBatch, row, col}); |
| 545 | }; |
| 546 | |
| 547 | bool converged = false; |
| 548 | |
| 549 | for (int i = 0; i < kMaxInterpolationSteps; i++) |
| 550 | { |
| 551 | // clang-format off |
| 552 | dD[0] = (dogVal(l + 0, r + 0, c + 1) - dogVal(l + 0, r + 0, c - 1)) * kDScale1; |
| 553 | dD[1] = (dogVal(l + 0, r + 1, c + 0) - dogVal(l + 0, r - 1, c + 0)) * kDScale1; |
| 554 | dD[2] = (dogVal(l + 1, r + 0, c + 0) - dogVal(l - 1, r + 0, c + 0)) * kDScale1; |
| 555 | |
| 556 | cv = dogVal(l, r, c); |
| 557 | |
| 558 | H[0][0] = (dogVal(l + 0, r + 0, c + 1) + dogVal(l + 0, r + 0, c - 1) - 2 * cv) * kDScale2; |
| 559 | H[1][1] = (dogVal(l + 0, r + 1, c + 0) + dogVal(l + 0, r - 1, c + 0) - 2 * cv) * kDScale2; |
| 560 | H[2][2] = (dogVal(l + 1, r + 0, c + 0) + dogVal(l - 1, r + 0, c + 0) - 2 * cv) * kDScale2; |
| 561 | |
| 562 | H[0][1] = H[1][0] = (dogVal(l + 0, r + 1, c + 1) - dogVal(l + 0, r + 1, c - 1) - |
| 563 | dogVal(l + 0, r - 1, c + 1) + dogVal(l + 0, r - 1, c - 1)) * kDScaleC; |
| 564 | H[0][2] = H[2][0] = (dogVal(l + 1, r + 0, c + 1) - dogVal(l + 1, r + 0, c - 1) - |
| 565 | dogVal(l - 1, r + 0, c + 1) + dogVal(l - 1, r + 0, c - 1)) * kDScaleC; |
| 566 | H[1][2] = H[2][1] = (dogVal(l + 1, r + 1, c + 0) - dogVal(l + 1, r - 1, c + 0) - |
| 567 | dogVal(l - 1, r + 1, c + 0) + dogVal(l - 1, r - 1, c + 0)) * kDScaleC; |
| 568 | // clang-format on |
| 569 | |
| 570 | sol = dD; |
| 571 | |
| 572 | if (!cuda::math::solve_inplace(H, sol)) |
| 573 | { |
| 574 | return; |
| 575 | } |
| 576 | |
| 577 | sol = -sol; |
| 578 | |
| 579 | if (std::abs(sol[2]) < 0.5f && std::abs(sol[1]) < 0.5f && std::abs(sol[0]) < 0.5f) |
| 580 | { |
| 581 | converged = true; |
| 582 | break; |
no test coverage detected