Skip to content

Commit 730baeb

Browse files
lusorisclaude
andcommitted
feature/adm/cuda: clamp the scale 0 contrast masking neighbours at the border
adm_cm_line_kernel() reads the neighbours of a sample at positions abs(x - 1), x and x + 1, and subtracts max(0, 2 * (x - w) + 1) from the last one (rows likewise). The correction tests the sample's own position against the band, not the neighbour's, so it never applies: x is inside the band. For a band of 14 samples or fewer in a dimension, a frame of 17 to 28 pixels, the contrast masking runs to the band's last column and row, and the neighbour at x + 1 is column w (row h), outside the band. The CPU replicates the last sample: {w - 2, w - 1, w - 1}. Clamp both neighbours with min(position, n - 1). The rows of a thread that lie below the band are masked off and unchanged, but they no longer read outside it. Add test_cuda_adm_cm_scale0_border, which compares integer_adm_scale0 of adm_cuda with the scalar CPU adm on frames of 17 to 28 pixels at 8 and 10 bits. Applies on top of the row rounding change: that one removes the rounding error that otherwise hides this one. Co-Authored-By: Claude Opus 5.5 <[email protected]>
1 parent 32884ee commit 730baeb

3 files changed

Lines changed: 278 additions & 4 deletions

File tree

‎libvmaf/src/feature/cuda/integer_adm/adm_cm.cu‎

Lines changed: 6 additions & 4 deletions
Original file line numberDiff line numberDiff line change
@@ -144,8 +144,11 @@ __device__ __forceinline__ void adm_cm_line_kernel(AdmBufferCuda buf, int h, int
144144
if (y < end_row && x < end_col)
145145
{
146146
int pos_x[3] = {x - 1, x, x + 1};
147-
pos_x[0] = abs(pos_x[0]);
148-
pos_x[2] = pos_x[2] - max(0, 2*(x - w)+1);
147+
// the left border mirrors {1, 0, 1} and the right border replicates
148+
// {w - 2, w - 1, w - 1}, as the CPU's ADM_CM_THRESH_S_* macros do. Clamp the
149+
// neighbour's position: x itself is inside the band.
150+
pos_x[0] = min(abs(pos_x[0]), w - 1);
151+
pos_x[2] = min(pos_x[2], w - 1);
149152

150153
#pragma unroll
151154
for (int theta = 0; theta < 3; ++theta)
@@ -165,8 +168,7 @@ __device__ __forceinline__ void adm_cm_line_kernel(AdmBufferCuda buf, int h, int
165168
for (int row = 0; row < total_rows;++row)
166169
{
167170
int pos_y = y - 1 + row;
168-
pos_y = abs(pos_y);
169-
pos_y = pos_y - max(0, 2*(y - h)+1);
171+
pos_y = min(abs(pos_y), h - 1);
170172

171173
int16_t src = angles[theta][pos_y * src_stride + x];
172174
int16_t *flt_ptr = flt_angles[theta] + pos_y*src_stride;

‎libvmaf/test/meson.build‎

Lines changed: 9 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -204,9 +204,18 @@ test_cuda_adm_cm_row_rounding = executable('test_cuda_adm_cm_row_rounding',
204204
c_args: ['-DHAVE_CUDA=1']
205205
)
206206

207+
test_cuda_adm_cm_scale0_border = executable('test_cuda_adm_cm_scale0_border',
208+
['test.c', 'test_cuda_adm_cm_scale0_border.c'],
209+
include_directories : [libvmaf_inc, test_inc],
210+
link_with : get_option('default_library') == 'both' ? libvmaf.get_static_lib() : libvmaf,
211+
dependencies: [math_lib, cuda_dependency],
212+
c_args: ['-DHAVE_CUDA=1']
213+
)
214+
207215
test('test_ring_buffer', test_ring_buffer)
208216
test('test_cuda_pic_preallocation', test_cuda_pic_preallocation)
209217
test('test_cuda_adm_cm_row_rounding', test_cuda_adm_cm_row_rounding, timeout : 300)
218+
test('test_cuda_adm_cm_scale0_border', test_cuda_adm_cm_scale0_border, timeout : 300)
210219
endif
211220

212221
test_pic_preallocation = executable('test_pic_preallocation',
Lines changed: 263 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,263 @@
1+
/**
2+
*
3+
* Copyright 2016-2023 Netflix, Inc.
4+
*
5+
* Licensed under the BSD+Patent License (the "License");
6+
* you may not use this file except in compliance with the License.
7+
* You may obtain a copy of the License at
8+
*
9+
* https://opensource.org/licenses/BSDplusPatent
10+
*
11+
* Unless required by applicable law or agreed to in writing, software
12+
* distributed under the License is distributed on an "AS IS" BASIS,
13+
* WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
14+
* See the License for the specific language governing permissions and
15+
* limitations under the License.
16+
*
17+
*/
18+
19+
/* integer_adm_scale0 of adm_cuda against the scalar CPU adm on frames 17 to 28
20+
* pixels wide or high, where the scale 0 band is 14 samples or fewer and the
21+
* contrast masking reads its right and bottom neighbours at the border: they
22+
* are replicated, {n - 2, n - 1, n - 1}. Needs a CUDA device. */
23+
24+
#include <math.h>
25+
#include <stdint.h>
26+
#include <stdio.h>
27+
#include <stdlib.h>
28+
#include <string.h>
29+
30+
#include "test.h"
31+
32+
#include "libvmaf/feature.h"
33+
#include "libvmaf/libvmaf.h"
34+
#include "libvmaf/libvmaf_cuda.h"
35+
#include "libvmaf/picture.h"
36+
37+
#define MAX_FRAMES 4
38+
#define NUM_SCORES 5
39+
40+
enum content {
41+
CONTENT_RANDOM,
42+
CONTENT_STRUCTURED,
43+
};
44+
45+
static const char *const score_name[NUM_SCORES] = {
46+
"VMAF_integer_feature_adm2_score",
47+
"integer_adm_scale0",
48+
"integer_adm_scale1",
49+
"integer_adm_scale2",
50+
"integer_adm_scale3",
51+
};
52+
53+
static uint64_t rng_state;
54+
55+
static uint32_t rng_next(void) {
56+
rng_state = rng_state * 6364136223846793005ULL + 1442695040888963407ULL;
57+
return (uint32_t)(rng_state >> 33);
58+
}
59+
60+
static int clip_sample(int v, int max) {
61+
return v < 0 ? 0 : (v > max ? max : v);
62+
}
63+
64+
/* Plane of the reference picture. Random: independent samples. Structured: a
65+
* ramp with a 3x3 checkerboard and a mirrored lower right quadrant, whose
66+
* smooth areas leave near-zero contrast masking accumulators. */
67+
static void fill_reference(int *ref, unsigned w, unsigned h, int mx,
68+
enum content content) {
69+
for (unsigned y = 0; y < h; y++) {
70+
for (unsigned x = 0; x < w; x++) {
71+
int v;
72+
if (content == CONTENT_RANDOM) {
73+
v = (int)(rng_next() % (unsigned)(mx + 1));
74+
} else {
75+
const int ramp = ((int)(x * (unsigned)mx / (w > 1 ? w - 1 : 1)) +
76+
(int)(y * (unsigned)mx / (h > 1 ? h - 1 : 1))) /
77+
2;
78+
const int check = (((x / 3) + (y / 3)) & 1) ? mx / 8 : 0;
79+
v = ramp / 2 + check + (int)(rng_next() % (unsigned)(mx / 32 + 1));
80+
if (x > w / 2 && y > h / 3)
81+
v = mx - v;
82+
}
83+
ref[y * w + x] = clip_sample(v, mx);
84+
}
85+
}
86+
}
87+
88+
/* Distorted plane: random noise on the reference, or its 3x3 box blur plus
89+
* a little noise. */
90+
static void fill_distorted(const int *ref, int *dis, unsigned w, unsigned h,
91+
int mx, enum content content) {
92+
for (unsigned y = 0; y < h; y++) {
93+
for (unsigned x = 0; x < w; x++) {
94+
int v;
95+
if (content == CONTENT_RANDOM) {
96+
v = ref[y * w + x] + (int)(rng_next() % (unsigned)(mx / 4 + 1)) -
97+
mx / 8;
98+
} else {
99+
int sum = 0, n = 0;
100+
for (int dy = -1; dy <= 1; dy++) {
101+
for (int dx = -1; dx <= 1; dx++) {
102+
const int yy = (int)y + dy, xx = (int)x + dx;
103+
if (yy < 0 || xx < 0 || yy >= (int)h || xx >= (int)w)
104+
continue;
105+
sum += ref[yy * w + xx];
106+
n++;
107+
}
108+
}
109+
v = sum / n + (int)(rng_next() % (unsigned)(mx / 64 + 1));
110+
}
111+
dis[y * w + x] = clip_sample(v, mx);
112+
}
113+
}
114+
}
115+
116+
static void store_plane(VmafPicture *pic, unsigned plane, const int *src) {
117+
const unsigned w = pic->w[plane], h = pic->h[plane];
118+
for (unsigned y = 0; y < h; y++) {
119+
uint8_t *row = (uint8_t *)pic->data[plane] + y * pic->stride[plane];
120+
for (unsigned x = 0; x < w; x++) {
121+
if (pic->bpc == 8)
122+
row[x] = (uint8_t)src[y * w + x];
123+
else
124+
((uint16_t *)row)[x] = (uint16_t)src[y * w + x];
125+
}
126+
}
127+
}
128+
129+
static int make_pictures(VmafPicture *ref, VmafPicture *dis, unsigned w,
130+
unsigned h, unsigned bpc, enum content content) {
131+
int err = vmaf_picture_alloc(ref, VMAF_PIX_FMT_YUV420P, bpc, w, h);
132+
err |= vmaf_picture_alloc(dis, VMAF_PIX_FMT_YUV420P, bpc, w, h);
133+
if (err)
134+
return err;
135+
const int mx = (1 << bpc) - 1;
136+
for (unsigned p = 0; p < 3; p++) {
137+
const unsigned n = ref->w[p] * ref->h[p];
138+
int *r = malloc(n * sizeof(*r));
139+
int *d = malloc(n * sizeof(*d));
140+
if (!r || !d) {
141+
free(r);
142+
free(d);
143+
return -1;
144+
}
145+
fill_reference(r, ref->w[p], ref->h[p], mx, content);
146+
fill_distorted(r, d, ref->w[p], ref->h[p], mx, content);
147+
store_plane(ref, p, r);
148+
store_plane(dis, p, d);
149+
free(r);
150+
free(d);
151+
}
152+
return 0;
153+
}
154+
155+
/* Scores of `frames` frames, from the CUDA extractor or from the scalar CPU
156+
* extractor. */
157+
static int run_adm(int use_cuda, unsigned w, unsigned h, unsigned bpc,
158+
enum content content, unsigned seed, unsigned frames,
159+
double scores[MAX_FRAMES][NUM_SCORES]) {
160+
VmafConfiguration cfg = {0};
161+
if (!use_cuda)
162+
cfg.cpumask = ~(uint64_t)0; /* scalar C only */
163+
164+
VmafContext *vmaf;
165+
if (vmaf_init(&vmaf, cfg))
166+
return -1;
167+
168+
int err = 0;
169+
if (use_cuda) {
170+
VmafCudaState *cu_state;
171+
VmafCudaConfiguration cu_cfg = {0};
172+
err = vmaf_cuda_state_init(&cu_state, cu_cfg);
173+
if (!err)
174+
err = vmaf_cuda_import_state(vmaf, cu_state);
175+
}
176+
if (!err)
177+
err = vmaf_use_feature(vmaf, use_cuda ? "adm_cuda" : "adm", NULL);
178+
179+
rng_state = (uint64_t)seed * 7919 + w * 131 + h;
180+
for (unsigned i = 0; !err && i < frames; i++) {
181+
VmafPicture ref, dis;
182+
err = make_pictures(&ref, &dis, w, h, bpc, content);
183+
if (!err)
184+
err = vmaf_read_pictures(vmaf, &ref, &dis, i);
185+
}
186+
if (!err)
187+
err = vmaf_read_pictures(vmaf, NULL, NULL, 0);
188+
for (unsigned i = 0; !err && i < frames; i++) {
189+
for (unsigned s = 0; s < NUM_SCORES; s++)
190+
err |= vmaf_feature_score_at_index(vmaf, score_name[s], &scores[i][s], i);
191+
}
192+
err |= vmaf_close(vmaf);
193+
return err;
194+
}
195+
196+
static int have_cuda_device(void) {
197+
VmafContext *vmaf;
198+
VmafConfiguration cfg = {0};
199+
if (vmaf_init(&vmaf, cfg))
200+
return 0;
201+
VmafCudaState *cu_state;
202+
VmafCudaConfiguration cu_cfg = {0};
203+
const int ok = !vmaf_cuda_state_init(&cu_state, cu_cfg);
204+
vmaf_close(vmaf);
205+
return ok;
206+
}
207+
208+
/* Largest absolute difference between the CUDA and the CPU score over the
209+
* frames, for the scores first..last. */
210+
static int max_difference(unsigned w, unsigned h, unsigned bpc,
211+
enum content content, unsigned seed, unsigned frames,
212+
unsigned first, unsigned last, double *worst) {
213+
double cuda[MAX_FRAMES][NUM_SCORES], cpu[MAX_FRAMES][NUM_SCORES];
214+
int err = run_adm(1, w, h, bpc, content, seed, frames, cuda);
215+
err |= run_adm(0, w, h, bpc, content, seed, frames, cpu);
216+
if (err)
217+
return err;
218+
*worst = 0.;
219+
for (unsigned i = 0; i < frames; i++) {
220+
for (unsigned s = first; s <= last; s++) {
221+
const double d = fabs(cuda[i][s] - cpu[i][s]);
222+
if (!(d <= *worst))
223+
*worst = d; /* also records NaN */
224+
}
225+
}
226+
return 0;
227+
}
228+
229+
#define TOLERANCE 1e-9
230+
231+
static const struct {
232+
unsigned w, h;
233+
} sizes[] = {{17, 17}, {20, 20}, {24, 24}, {28, 28}, {17, 64}, {64, 17}};
234+
235+
static char *test_scale_0_border_matches_cpu(void)
236+
{
237+
for (unsigned i = 0; i < sizeof(sizes) / sizeof(sizes[0]); i++) {
238+
for (unsigned bpc = 8; bpc <= 10; bpc += 2) {
239+
for (int content = CONTENT_RANDOM; content <= CONTENT_STRUCTURED; content++) {
240+
double worst;
241+
const int err = max_difference(sizes[i].w, sizes[i].h, bpc, content, 1, 2, 1, 1,
242+
&worst);
243+
mu_assert("problem during the CUDA / CPU comparison", !err);
244+
if (worst > TOLERANCE)
245+
fprintf(stderr, "%ux%u %u-bit content %d: |cuda - cpu| = %g\n", sizes[i].w,
246+
sizes[i].h, bpc, content, worst);
247+
mu_assert("integer_adm_scale0 of adm_cuda differs from the CPU", worst <= TOLERANCE);
248+
}
249+
}
250+
}
251+
return NULL;
252+
}
253+
254+
char *run_tests(void)
255+
{
256+
if (!have_cuda_device()) {
257+
fprintf(stderr, "no CUDA device: skipped\n");
258+
return NULL;
259+
}
260+
mu_run_test(test_scale_0_border_matches_cpu);
261+
return NULL;
262+
}
263+

0 commit comments

Comments
 (0)