48 |
|
|
49 |
typedef struct { |
typedef struct { |
50 |
|
|
|
int pixels; /* width*height */ |
|
51 |
uint64_t mse_sum; /* for avrg psnrhvsm */ |
uint64_t mse_sum; /* for avrg psnrhvsm */ |
52 |
long frame_cnt; |
long frame_cnt; |
53 |
|
|
75 |
0.019290f, 0.011815f, 0.011080f, 0.010412f, 0.007972f, 0.010000f, 0.009426f, 0.010203f |
0.019290f, 0.011815f, 0.011080f, 0.010412f, 0.007972f, 0.010000f, 0.009426f, 0.010203f |
76 |
}; |
}; |
77 |
|
|
78 |
#if 1 /* Floating-point implementation */ |
#if 0 /* Floating-point implementation */ |
79 |
|
|
80 |
static uint32_t Calc_MSE_H(int16_t *DCT_A, int16_t *DCT_B, uint8_t *IMG_A, uint8_t *IMG_B, int stride) |
static uint32_t Calc_MSE_H(int16_t *DCT_A, int16_t *DCT_B, uint8_t *IMG_A, uint8_t *IMG_B, int stride) |
81 |
{ |
{ |
84 |
uint32_t Local[8] = {0, 0, 0, 0, 0, 0, 0, 0}; |
uint32_t Local[8] = {0, 0, 0, 0, 0, 0, 0, 0}; |
85 |
uint32_t Local_Square[8] = {0, 0, 0, 0, 0, 0, 0, 0}; |
uint32_t Local_Square[8] = {0, 0, 0, 0, 0, 0, 0, 0}; |
86 |
float MASK_A = 0.f, MASK_B = 0.f; |
float MASK_A = 0.f, MASK_B = 0.f; |
87 |
|
float Mult1 = 1.f, Mult2 = 1.f; |
88 |
uint32_t MSE_H = 0; |
uint32_t MSE_H = 0; |
89 |
|
|
90 |
/* Step 1: Calculate CSF weighted energy of DCT coefficients */ |
/* Step 1: Calculate CSF weighted energy of DCT coefficients */ |
98 |
/* Step 2: Determine local variances compared to entire block variance */ |
/* Step 2: Determine local variances compared to entire block variance */ |
99 |
for (y = 0; y < 2; y++) { |
for (y = 0; y < 2; y++) { |
100 |
for (x = 0; x < 2; x++) { |
for (x = 0; x < 2; x++) { |
101 |
|
for (j = 0; j < 4; j++) { |
102 |
for (i = 0; i < 4; i++) { |
for (i = 0; i < 4; i++) { |
103 |
uint8_t A = IMG_A[y*4*stride + 4*x + i]; |
uint8_t A = IMG_A[(y*4+j)*stride + 4*x + i]; |
104 |
uint8_t B = IMG_B[y*4*stride + 4*x + i]; |
uint8_t B = IMG_B[(y*4+j)*stride + 4*x + i]; |
105 |
|
|
106 |
Local[y*2 + x] += A; |
Local[y*2 + x] += A; |
107 |
Local[y*2 + x + 4] += B; |
Local[y*2 + x + 4] += B; |
110 |
} |
} |
111 |
} |
} |
112 |
} |
} |
113 |
|
} |
114 |
|
|
115 |
Global_A = Local[0] + Local[1] + Local[2] + Local[3]; |
Global_A = Local[0] + Local[1] + Local[2] + Local[3]; |
116 |
Global_B = Local[4] + Local[5] + Local[6] + Local[7]; |
Global_B = Local[4] + Local[5] + Local[6] + Local[7]; |
125 |
Global_B = (Local_Square[4]<<6) - Global_B*Global_B; /* 64*Var(D) */ |
Global_B = (Local_Square[4]<<6) - Global_B*Global_B; /* 64*Var(D) */ |
126 |
|
|
127 |
/* Step 3: Calculate contrast masking threshold */ |
/* Step 3: Calculate contrast masking threshold */ |
128 |
MASK_A = (float)sqrt(MASK_A*((float)(Local[0]+Local[1]+Local[2]+Local[3])/((float)Global_A/4.f)))/32.f; |
if (Global_A) |
129 |
MASK_B = (float)sqrt(MASK_B*((float)(Local[4]+Local[5]+Local[6]+Local[7])/((float)Global_B/4.f)))/32.f; |
Mult1 = (float)(Local[0]+Local[1]+Local[2]+Local[3])/((float)(Global_A)/4.f); |
130 |
|
|
131 |
|
if (Global_B) |
132 |
|
Mult2 = (float)(Local[4]+Local[5]+Local[6]+Local[7])/((float)(Global_B)/4.f); |
133 |
|
|
134 |
|
MASK_A = (float)sqrt(MASK_A * Mult1) / 32.f; |
135 |
|
MASK_B = (float)sqrt(MASK_B * Mult2) / 32.f; |
136 |
|
|
137 |
if (MASK_B > MASK_A) MASK_A = MASK_B; /* MAX(MASK_A, MASK_B) */ |
if (MASK_B > MASK_A) MASK_A = MASK_B; /* MAX(MASK_A, MASK_B) */ |
138 |
|
|
190 |
366, 286, 277, 269, 235, 264, 256, 266 |
366, 286, 277, 269, 235, 264, 256, 266 |
191 |
}; |
}; |
192 |
|
|
193 |
|
static __inline uint32_t isqrt(unsigned long n) |
194 |
|
{ |
195 |
|
uint32_t c = 0x8000; |
196 |
|
uint32_t g = 0x8000; |
197 |
|
|
198 |
|
for(;;) { |
199 |
|
if(g*g > n) |
200 |
|
g ^= c; |
201 |
|
c >>= 1; |
202 |
|
if(c == 0) |
203 |
|
return g; |
204 |
|
g |= c; |
205 |
|
} |
206 |
|
} |
207 |
|
|
208 |
static uint32_t Calc_MSE_H(int16_t *DCT_A, int16_t *DCT_B, uint8_t *IMG_A, uint8_t *IMG_B, int stride) |
static uint32_t Calc_MSE_H(int16_t *DCT_A, int16_t *DCT_B, uint8_t *IMG_A, uint8_t *IMG_B, int stride) |
209 |
{ |
{ |
210 |
int x, y, i, j; |
int x, y, i, j; |
228 |
/* Step 2: Determine local variances compared to entire block variance */ |
/* Step 2: Determine local variances compared to entire block variance */ |
229 |
for (y = 0; y < 2; y++) { |
for (y = 0; y < 2; y++) { |
230 |
for (x = 0; x < 2; x++) { |
for (x = 0; x < 2; x++) { |
231 |
|
for (j = 0; j < 4; j++) { |
232 |
for (i = 0; i < 4; i++) { |
for (i = 0; i < 4; i++) { |
233 |
uint8_t A = IMG_A[y*4*stride + 4*x + i]; |
uint8_t A = IMG_A[(y*4+j)*stride + 4*x + i]; |
234 |
uint8_t B = IMG_B[y*4*stride + 4*x + i]; |
uint8_t B = IMG_B[(y*4+j)*stride + 4*x + i]; |
235 |
|
|
236 |
Local[y*2 + x] += A; |
Local[y*2 + x] += A; |
237 |
Local[y*2 + x + 4] += B; |
Local[y*2 + x + 4] += B; |
240 |
} |
} |
241 |
} |
} |
242 |
} |
} |
243 |
|
} |
244 |
|
|
245 |
Global_A = Local[0] + Local[1] + Local[2] + Local[3]; |
Global_A = Local[0] + Local[1] + Local[2] + Local[3]; |
246 |
Global_B = Local[4] + Local[5] + Local[6] + Local[7]; |
Global_B = Local[4] + Local[5] + Local[6] + Local[7]; |
256 |
|
|
257 |
/* Step 3: Calculate contrast masking threshold */ |
/* Step 3: Calculate contrast masking threshold */ |
258 |
{ |
{ |
259 |
float MASK_A, MASK_B; /* TODO */ |
uint32_t MASK_A, MASK_B; |
260 |
|
uint32_t Mult1 = 64, Mult2 = 64; |
261 |
|
|
262 |
|
if (Global_A) |
263 |
|
Mult1 = ((Local[0]+Local[1]+Local[2]+Local[3])<<8) / Global_A; |
264 |
|
|
265 |
MASK_A = (float)((double)(Sum_A + 4) / 8.); |
if (Global_B) |
266 |
MASK_B = (float)((double)(Sum_B + 4) / 8.); |
Mult2 = ((Local[4]+Local[5]+Local[6]+Local[7])<<8) / Global_B; |
267 |
|
|
268 |
MASK_A = sqrt(MASK_A*((float)(Local[0]+Local[1]+Local[2]+Local[3])/((float)Global_A/4.f)))/32.f; |
MASK_A = isqrt(2*Sum_A*Mult1) + 16; |
269 |
MASK_B = sqrt(MASK_B*((float)(Local[4]+Local[5]+Local[6]+Local[7])/((float)Global_B/4.f)))/32.f; |
MASK_B = isqrt(2*Sum_B*Mult2) + 16; |
270 |
|
|
271 |
if (MASK_B > MASK_A) /* MAX(MASK_A, MASK_B) */ |
if (MASK_B > MASK_A) /* MAX(MASK_A, MASK_B) */ |
272 |
MASK = (uint32_t) (MASK_B * 1024.f + 0.5f); |
MASK = (uint32_t) MASK_B; |
273 |
else |
else |
274 |
MASK = (uint32_t) (MASK_A * 1024.f + 0.5f); |
MASK = (uint32_t) MASK_A; |
275 |
} |
} |
276 |
|
|
277 |
/* Step 4: Calculate MSE of DCT coeffs reduced by masking effect */ |
/* Step 4: Calculate MSE of DCT coeffs reduced by masking effect */ |
301 |
static void psnrhvsm_after(xvid_plg_data_t *data, psnrhvsm_data_t *psnrhvsm) |
static void psnrhvsm_after(xvid_plg_data_t *data, psnrhvsm_data_t *psnrhvsm) |
302 |
{ |
{ |
303 |
DECLARE_ALIGNED_MATRIX(DCT, 2, 64, int16_t, CACHE_LINE); |
DECLARE_ALIGNED_MATRIX(DCT, 2, 64, int16_t, CACHE_LINE); |
304 |
int x, y, stride = data->original.stride[0]; |
uint32_t x, y, stride = data->original.stride[0]; |
305 |
int16_t *DCT_A = &DCT[0], *DCT_B = &DCT[64]; |
int16_t *DCT_A = &DCT[0], *DCT_B = &DCT[64]; |
306 |
uint8_t *IMG_A = (uint8_t *) data->original.plane[0]; |
uint8_t *IMG_A = (uint8_t *) data->original.plane[0]; |
307 |
uint8_t *IMG_B = (uint8_t *) data->current.plane[0]; |
uint8_t *IMG_B = (uint8_t *) data->current.plane[0]; |
308 |
uint32_t MSE_H = 0; |
uint64_t MSE_H = 0; |
|
|
|
|
psnrhvsm->pixels = data->width * data->height; |
|
309 |
|
|
310 |
for (y = 0; y < data->height; y += 8) { |
for (y = 0; y < data->height; y += 8) { |
311 |
for (x = 0; x < data->width; x += 8) { |
for (x = 0; x < data->width; x += 8) { |
328 |
} |
} |
329 |
} |
} |
330 |
|
|
331 |
psnrhvsm->mse_sum += MSE_H; |
x = 4*MSE_H / (data->width * data->height); |
332 |
|
psnrhvsm->mse_sum += x; |
333 |
psnrhvsm->frame_cnt++; |
psnrhvsm->frame_cnt++; |
334 |
|
|
335 |
printf(" psnrhvsm: %2.2f\n", sse_to_PSNR(MSE_H, 256*psnrhvsm->pixels)); |
printf(" psnrhvsm: %2.2f\n", sse_to_PSNR(x, 1024)); |
336 |
} |
} |
337 |
|
|
338 |
static int psnrhvsm_create(xvid_plg_create_t *create, void **handle) |
static int psnrhvsm_create(xvid_plg_create_t *create, void **handle) |
343 |
psnrhvsm->mse_sum = 0; |
psnrhvsm->mse_sum = 0; |
344 |
psnrhvsm->frame_cnt = 0; |
psnrhvsm->frame_cnt = 0; |
345 |
|
|
|
psnrhvsm->pixels = 0; |
|
|
|
|
346 |
*(handle) = (void*) psnrhvsm; |
*(handle) = (void*) psnrhvsm; |
347 |
return 0; |
return 0; |
348 |
} |
} |
371 |
MSE_H = (uint32_t) (psnrhvsm->mse_sum / psnrhvsm->frame_cnt); |
MSE_H = (uint32_t) (psnrhvsm->mse_sum / psnrhvsm->frame_cnt); |
372 |
|
|
373 |
emms(); |
emms(); |
374 |
printf("Average psnrhvsm: %2.2f\n", sse_to_PSNR(MSE_H, 256*psnrhvsm->pixels)); |
printf("Average psnrhvsm: %2.2f\n", sse_to_PSNR(MSE_H, 1024)); |
375 |
free(psnrhvsm); |
free(psnrhvsm); |
376 |
} |
} |
377 |
} |
} |