potok/opencv
1
0
Fork 0
mirror of https://github.com/opencv/opencv.git synced 2026-07-21 19:33:03 +04:00
Code Packages Releases Wiki Activity

Merge pull request #29194 from Enilrats:opt-emd-performance

Browse Source 
imgproc: optimize EMD (Earth Mover's Distance) solver performance #29194

### imgproc: optimize EMD solver performance using O(V) spanning tree traversal

---
This PR significantly optimizes the performance of the Earth Mover's Distance (`cv::EMD`) solver in `modules/imgproc/src/emd_new.cpp`. 

Specifically, it refactors the dual-variable calculation in `EMDSolver::findBasicVars()` from a naive $O((N+M)^2)$ linked-list scanning approach to an optimized $O(N+M)$ BFS tree traversal utilizing existing adjacency lists, leading to a massive speedup.

---

#### Technical Details & Core Bottleneck Fixed

1. **Algorithmic Complexity Reduction in `findBasicVars()`:**
   - **Before**: The original implementation solved the dual variables $u_i$ and $v_j$ by traversing the entire unmarked rows (`u0_head`) or columns (`v0_head`) linked-lists and invoking `getIsX(i, j)` inside nested loops to find connected basic variables. This resulted in an $O((N+M)^2)$ complexity per simplex iteration. For a scale of $2000 \times 2000$, this performed up to $16,000,000$ operations per iteration.
   - **After**: Since `EMDSolver` already maintains the adjacency lists of the basic variables tree (`rows_x` and `cols_x`), we can traverse the spanning tree in linear time. This PR implements a dual-queue BFS tree traversal. The complexity per iteration is drastically reduced to $O(N+M)$, performing at most $4000$ operations per iteration.

2. **Cache Locality & Pointer-Chasing Elimination:**
   - Replaced pointer-chasing on dynamically-allocated linked lists with contiguous, stack-allocated array queues (`cv::AutoBuffer`), significantly improving CPU L1/L2 cache hit rates and enabling hardware prefetching.


---

#### Performance Benchmarks

Below is the benchmark comparison evaluated on a standard CPU.
#### Test 1: dims = 64
| Scale ($N, M$) | Original EMD (ms) | Optimized EMD (This PR) | Speedup |
| :--- | :--- | :--- | :--- |
| **100** | 5.473 | 3.419 | **1.60x** |
| **500** | 381.871 | 244.355 | **1.56x** |
| **1000** | 1893.053 | 1369.016 | **1.38x** |
| **2000** | 11387.792 | 8331.221 | **1.37x** |

#### Test 2: dims = 3
| Scale ($N, M$) | Original EMD (ms) | Optimized EMD (This PR) | Speedup |
| :--- | :--- | :--- | :--- |
| **100** | 4.433 | 3.042 | **1.46x** |
| **500** | 365.762 | 259.735 | **1.41x** |
| **1000** | 1989.400 | 1421.952 | **1.40x** |
| **2000** | 12731.836 | 7952.210 | **1.60x** |


*(Note: The exact performance figures may vary slightly depending on the compiler and test machine.)*

---

#### Verification
- All existing tests in `opencv_test_imgproc` (including EMD tests) pass successfully. No regressions were introduced.

---

### Pull Request Readiness Checklist

See details at https://github.com/opencv/opencv/wiki/How_to_contribute#making-a-good-pull-request

- [x] I agree to contribute to the project under Apache 2 License.
- [x] To the best of my knowledge, the proposed patch is not based on a code under GPL or another license that is incompatible with OpenCV
- [x] The PR is proposed to the proper branch
- [ ] There is a reference to the original bug report and related work
- [ ] There is accuracy test, performance test and test data in opencv_extra repository, if applicable
      Patch to opencv_extra has the same branch name.
- [ ] The feature is well documented and sample code can be built with the project CMake
This commit is contained in:
Ding zhehao
2026-06-10 23:23:35 +08:00
committed by GitHub
parent 1780a86075
commit ffa38e1b74
2 changed files with 88 additions and 92 deletions
Download Patch File Download Diff File Expand all files Collapse all files
modules/imgproc/perf/perf_emd.cpp
+37
View File
@@ -0,0 +1,37 @@
// This file is part of OpenCV project.
// It is subject to the license terms in the LICENSE file found in the top-level directory
// of this distribution and at http://opencv.org/license.html
#include "perf_precomp.hpp"
#include <cmath>
namespace opencv_test { namespace {
typedef tuple<int, int> EMD_Size_Dim_t;
typedef perf::TestBaseWithParam<EMD_Size_Dim_t> EMD_Fixture;
PERF_TEST_P(EMD_Fixture, L1_Distance, testing::Combine(
testing::Values(100, 500, 1000),
testing::Values(3, 64)
))
{
int size = get<0>(GetParam());
int dims = get<1>(GetParam());
Mat sign1(size, dims + 1, CV_32FC1);
Mat sign2(size, dims + 1, CV_32FC1);
theRNG().fill(sign1, RNG::UNIFORM, 0.1, 1.0);
theRNG().fill(sign2, RNG::UNIFORM, 0.1, 1.0);
declare.in(sign1, sign2);
TEST_CYCLE()
{
cv::EMD(sign1, sign2, cv::DIST_L1);
}
SANITY_CHECK_NOTHING();
}
}}
modules/imgproc/src/emd_new.cpp
+51 -92
View File
@@ -442,119 +442,78 @@ double EMDSolver::calcFlow(Mat* flow_) const
int EMDSolver::findBasicVars() const
{
int i, j;
int u_cfound, v_cfound;
Node1D u0_head, u1_head, *cur_u, *prev_u;
Node1D v0_head, v1_head, *cur_v, *prev_v;
bool found;
CV_Assert(u != 0 && v != 0);
/* initialize the rows list (u) and the columns list (v) */
u0_head.next = u;
for (i = 0; i < ssize; i++)
{
u[i].next = u + i + 1;
}
u[ssize - 1].next = 0;
u1_head.next = 0;
// 1. Initialize status flags using contiguous memory to eliminate pointer chasing
AutoBuffer<char> computed_buf(ssize + dsize);
char* row_computed = computed_buf.data();
char* col_computed = computed_buf.data() + ssize;
memset(row_computed, 0, ssize + dsize);
v0_head.next = ssize > 1 ? v + 1 : 0;
for (i = 1; i < dsize; i++)
{
v[i].next = v + i + 1;
}
v[dsize - 1].next = 0;
v1_head.next = 0;
// 2. Create BFS queues
AutoBuffer<int> queue_buf(ssize + dsize);
int* row_queue = queue_buf.data();
int* col_queue = queue_buf.data() + ssize;
/* there are ssize+dsize variables but only ssize+dsize-1 independent equations,
so set v[0]=0 */
int row_head = 0, row_tail = 0;
int col_head = 0, col_tail = 0;
// Initial condition: enqueue column 0 as the root node and set its value to 0
v[0].val = 0;
v1_head.next = v;
v1_head.next->next = 0;
col_computed[0] = true;
col_queue[col_tail++] = 0;
/* loop until all variables are found */
u_cfound = v_cfound = 0;
while (u_cfound < ssize || v_cfound < dsize)
int u_cfound = 0;
int v_cfound = 1;
// 3. Dual-queue interactive BFS traversal over the spanning tree (Time Complexity: O(N + M))
while (row_head < row_tail || col_head < col_tail)
{
found = false;
if (v_cfound < dsize)
// Process currently marked columns to update their connected rows
while (col_head < col_tail)
{
/* loop over all marked columns */
prev_v = &v1_head;
cur_v = v1_head.next;
found = found || (cur_v != 0);
for (; cur_v != 0; cur_v = cur_v->next)
{
float cur_v_val = cur_v->val;
int j = col_queue[col_head++];
float cur_v_val = v[j].val;
j = (int)(cur_v - v);
/* find the variables in column j */
prev_u = &u0_head;
for (cur_u = u0_head.next; cur_u != 0;)
// Use adjacency list cols_x to directly access rows connected to column j, avoiding full scans
for (Node2D* xp = cols_x[j]; xp != 0; xp = xp->next[1])
{
int i = xp->i;
if (!row_computed[i])
{
i = (int)(cur_u - u);
if (getIsX(i, j))
{
/* compute u[i] */
cur_u->val = getCost(i, j) - cur_v_val;
/* ...and add it to the marked list */
prev_u->next = cur_u->next;
cur_u->next = u1_head.next;
u1_head.next = cur_u;
cur_u = prev_u->next;
}
else
{
prev_u = cur_u;
cur_u = cur_u->next;
}
u[i].val = getCost(i, j) - cur_v_val;
row_computed[i] = true;
row_queue[row_tail++] = i; // Enqueue the newly resolved row
u_cfound++;
}
prev_v->next = cur_v->next;
v_cfound++;
}
}
if (u_cfound < ssize)
// Process currently marked rows to update their connected columns
while (row_head < row_tail)
{
/* loop over all marked rows */
prev_u = &u1_head;
cur_u = u1_head.next;
found = found || (cur_u != 0);
for (; cur_u != 0; cur_u = cur_u->next)
int i = row_queue[row_head++];
float cur_u_val = u[i].val;
// Use adjacency list rows_x to directly access columns connected to row i
for (Node2D* xp = rows_x[i]; xp != 0; xp = xp->next[0])
{
float cur_u_val = cur_u->val;
i = (int)(cur_u - u);
/* find the variables in rows i */
prev_v = &v0_head;
for (cur_v = v0_head.next; cur_v != 0;)
int j = xp->j;
if (!col_computed[j])
{
j = (int)(cur_v - v);
if (getIsX(i, j))
{
/* compute v[j] */
cur_v->val = getCost(i, j) - cur_u_val;
/* ...and add it to the marked list */
prev_v->next = cur_v->next;
cur_v->next = v1_head.next;
v1_head.next = cur_v;
cur_v = prev_v->next;
}
else
{
prev_v = cur_v;
cur_v = cur_v->next;
}
v[j].val = getCost(i, j) - cur_u_val;
col_computed[j] = true;
col_queue[col_tail++] = j; // Enqueue the newly resolved column
v_cfound++;
}
prev_u->next = cur_u->next;
u_cfound++;
}
}
if (!found)
return -1;
}
// If the number of traversed nodes is insufficient, the graph is disconnected and the spanning tree is incomplete
if (u_cfound < ssize || v_cfound < dsize)
return -1;
return 0;
}
@@ -1008,4 +967,4 @@ float cv::wrapperEMD(InputArray _sign1,
OutputArray _flow)
{
return EMD(_sign1, _sign2, distType, _cost, lowerBound.get(), _flow);
}
}