tests/tiny_psnr.c
Go to the documentation of this file.
00001 /*
00002  * Copyright (c) 2003 Michael Niedermayer <michaelni@gmx.at>
00003  *
00004  * This file is part of FFmpeg.
00005  *
00006  * FFmpeg is free software; you can redistribute it and/or
00007  * modify it under the terms of the GNU Lesser General Public
00008  * License as published by the Free Software Foundation; either
00009  * version 2.1 of the License, or (at your option) any later version.
00010  *
00011  * FFmpeg is distributed in the hope that it will be useful,
00012  * but WITHOUT ANY WARRANTY; without even the implied warranty of
00013  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
00014  * Lesser General Public License for more details.
00015  *
00016  * You should have received a copy of the GNU Lesser General Public
00017  * License along with FFmpeg; if not, write to the Free Software
00018  * Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA
00019  */
00020 
00021 #include <stdio.h>
00022 #include <stdlib.h>
00023 #include <string.h>
00024 #include <inttypes.h>
00025 #include <assert.h>
00026 
00027 #define FFMIN(a, b) ((a) > (b) ? (b) : (a))
00028 #define F 100
00029 #define SIZE 2048
00030 
00031 uint64_t exp16_table[21] = {
00032            65537,
00033            65538,
00034            65540,
00035            65544,
00036            65552,
00037            65568,
00038            65600,
00039            65664,
00040            65793,
00041            66050,
00042            66568,
00043            67616,
00044            69763,
00045            74262,
00046            84150,
00047           108051,
00048           178145,
00049           484249,
00050          3578144,
00051        195360063,
00052     582360139072LL,
00053 };
00054 
00055 #if 0
00056 // 16.16 fixpoint exp()
00057 static unsigned int exp16(unsigned int a){
00058     int i;
00059     int out= 1<<16;
00060 
00061     for(i=19;i>=0;i--){
00062         if(a&(1<<i))
00063             out= (out*exp16_table[i] + (1<<15))>>16;
00064     }
00065 
00066     return out;
00067 }
00068 #endif
00069 
00070 // 16.16 fixpoint log()
00071 static int64_t log16(uint64_t a)
00072 {
00073     int i;
00074     int out = 0;
00075 
00076     if (a < 1 << 16)
00077         return -log16((1LL << 32) / a);
00078     a <<= 16;
00079 
00080     for (i = 20; i >= 0; i--) {
00081         int64_t b = exp16_table[i];
00082         if (a < (b << 16))
00083             continue;
00084         out |= 1 << i;
00085         a    = ((a / b) << 16) + (((a % b) << 16) + b / 2) / b;
00086     }
00087     return out;
00088 }
00089 
00090 static uint64_t int_sqrt(uint64_t a)
00091 {
00092     uint64_t ret    = 0;
00093     uint64_t ret_sq = 0;
00094     int s;
00095 
00096     for (s = 31; s >= 0; s--) {
00097         uint64_t b = ret_sq + (1ULL << (s * 2)) + (ret << s) * 2;
00098         if (b <= a) {
00099             ret_sq = b;
00100             ret   += 1ULL << s;
00101         }
00102     }
00103     return ret;
00104 }
00105 
00106 int main(int argc, char *argv[])
00107 {
00108     int i, j;
00109     uint64_t sse = 0;
00110     uint64_t dev;
00111     FILE *f[2];
00112     uint8_t buf[2][SIZE];
00113     uint64_t psnr;
00114     int len        = argc < 4 ? 1 : atoi(argv[3]);
00115     int64_t max    = (1 << (8 * len)) - 1;
00116     int shift      = argc < 5 ? 0 : atoi(argv[4]);
00117     int skip_bytes = argc < 6 ? 0 : atoi(argv[5]);
00118     int size0      = 0;
00119     int size1      = 0;
00120     int maxdist    = 0;
00121 
00122     if (argc < 3) {
00123         printf("tiny_psnr <file1> <file2> [<elem size> [<shift> [<skip bytes>]]]\n");
00124         printf("WAV headers are skipped automatically.\n");
00125         return 1;
00126     }
00127 
00128     f[0] = fopen(argv[1], "rb");
00129     f[1] = fopen(argv[2], "rb");
00130     if (!f[0] || !f[1]) {
00131         fprintf(stderr, "Could not open input files.\n");
00132         return 1;
00133     }
00134 
00135     for (i = 0; i < 2; i++) {
00136         uint8_t *p = buf[i];
00137         if (fread(p, 1, 12, f[i]) != 12)
00138             return 1;
00139         if (!memcmp(p, "RIFF", 4) &&
00140             !memcmp(p + 8, "WAVE", 4)) {
00141             if (fread(p, 1, 8, f[i]) != 8)
00142                 return 1;
00143             while (memcmp(p, "data", 4)) {
00144                 int s = p[4] | p[5] << 8 | p[6] << 16 | p[7] << 24;
00145                 fseek(f[i], s, SEEK_CUR);
00146                 if (fread(p, 1, 8, f[i]) != 8)
00147                     return 1;
00148             }
00149         } else {
00150             fseek(f[i], -12, SEEK_CUR);
00151         }
00152     }
00153 
00154     fseek(f[shift < 0], abs(shift), SEEK_CUR);
00155 
00156     fseek(f[0], skip_bytes, SEEK_CUR);
00157     fseek(f[1], skip_bytes, SEEK_CUR);
00158 
00159     for (;;) {
00160         int s0 = fread(buf[0], 1, SIZE, f[0]);
00161         int s1 = fread(buf[1], 1, SIZE, f[1]);
00162 
00163         for (j = 0; j < FFMIN(s0, s1); j++) {
00164             int64_t a = buf[0][j];
00165             int64_t b = buf[1][j];
00166             int dist;
00167             if (len == 2) {
00168                 a = (int16_t)(a | (buf[0][++j] << 8));
00169                 b = (int16_t)(b | (buf[1][  j] << 8));
00170             }
00171             sse += (a - b) * (a - b);
00172             dist = abs(a - b);
00173             if (dist > maxdist)
00174                 maxdist = dist;
00175         }
00176         size0 += s0;
00177         size1 += s1;
00178         if (s0 + s1 <= 0)
00179             break;
00180     }
00181 
00182     i = FFMIN(size0, size1) / len;
00183     if (!i)
00184         i = 1;
00185     dev = int_sqrt(((sse / i) * F * F) + (((sse % i) * F * F) + i / 2) / i);
00186     if (sse)
00187         psnr = ((2 * log16(max << 16) + log16(i) - log16(sse)) *
00188                 284619LL * F + (1LL << 31)) / (1LL << 32);
00189     else
00190         psnr = 1000 * F - 1; // floating point free infinity :)
00191 
00192     printf("stddev:%5d.%02d PSNR:%3d.%02d MAXDIFF:%5d bytes:%9d/%9d\n",
00193            (int)(dev / F), (int)(dev % F),
00194            (int)(psnr / F), (int)(psnr % F),
00195            maxdist, size0, size1);
00196     return 0;
00197 }