Statistics
| Branch: | Revision:

ffmpeg / tests / tiny_psnr.c @ ce611a27

History | View | Annotate | Download (3.61 KB)

1
/*
2
 * Copyright (c) 2003 Michael Niedermayer <michaelni@gmx.at>
3
 *
4
 * This file is part of FFmpeg.
5
 *
6
 * FFmpeg is free software; you can redistribute it and/or
7
 * modify it under the terms of the GNU Lesser General Public
8
 * License as published by the Free Software Foundation; either
9
 * version 2.1 of the License, or (at your option) any later version.
10
 *
11
 * FFmpeg is distributed in the hope that it will be useful,
12
 * but WITHOUT ANY WARRANTY; without even the implied warranty of
13
 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
14
 * Lesser General Public License for more details.
15
 *
16
 * You should have received a copy of the GNU Lesser General Public
17
 * License along with FFmpeg; if not, write to the Free Software
18
 * Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA
19
 */
20

    
21
#include <stdio.h>
22
#include <stdlib.h>
23
#include <inttypes.h>
24
#include <assert.h>
25

    
26
#define F 100
27
#define SIZE 2048
28

    
29
uint64_t exp16_table[21]={
30
     65537,
31
     65538,
32
     65540,
33
     65544,
34
     65552,
35
     65568,
36
     65600,
37
     65664,
38
     65793,
39
     66050,
40
     66568,
41
     67616,
42
     69763,
43
     74262,
44
     84150,
45
    108051,
46
    178145,
47
    484249,
48
   3578144,
49
 195360063,
50
 582360139072LL,
51
};
52
#if 1
53
// 16.16 fixpoint exp()
54
static unsigned int exp16(unsigned int a){
55
    int i;
56
    int out= 1<<16;
57

    
58
    for(i=19;i>=0;i--){
59
        if(a&(1<<i))
60
            out= (out*exp16_table[i] + (1<<15))>>16;
61
    }
62

    
63
    return out;
64
}
65
// 16.16 fixpoint log()
66
static int64_t log16(uint64_t a){
67
    int i;
68
    int out=0;
69

    
70
    if(a < 1<<16)
71
        return -log16((1LL<<32) / a);
72
    a<<=16;
73

    
74
    for(i=20;i>=0;i--){
75
        int64_t b= exp16_table[i];
76
        if(a<(b<<16)) continue;
77
        out |= 1<<i;
78
        a = ((a/b)<<16) + (((a%b)<<16) + b/2)/b;
79
    }
80
    return out;
81
}
82

    
83
#endif
84
static uint64_t int_sqrt(uint64_t a)
85
{
86
    uint64_t ret=0;
87
    int s;
88
    uint64_t ret_sq=0;
89

    
90
    for(s=31; s>=0; s--){
91
        uint64_t b= ret_sq + (1ULL<<(s*2)) + (ret<<s)*2;
92
        if(b<=a){
93
            ret_sq=b;
94
            ret+= 1ULL<<s;
95
        }
96
    }
97
    return ret;
98
}
99

    
100
int main(int argc,char* argv[]){
101
    int i, j;
102
    uint64_t sse=0;
103
    uint64_t dev;
104
    FILE *f[2];
105
    uint8_t buf[2][SIZE];
106
    uint64_t psnr;
107
    int len= argc<4 ? 1 : atoi(argv[3]);
108
    int64_t max= (1<<(8*len))-1;
109
    int shift= argc<5 ? 0 : atoi(argv[4]);
110
    int skip_bytes = argc<6 ? 0 : atoi(argv[5]);
111

    
112
    if(argc<3){
113
        printf("tiny_psnr <file1> <file2> [<elem size> [<shift> [<skip bytes>]]]\n");
114
        printf("for wav files use the following:\n");
115
        printf("./tiny_psnr file1.wav file2.wav 2 0 44 to skip the header.\n");
116
        return -1;
117
    }
118

    
119
    f[0]= fopen(argv[1], "rb");
120
    f[1]= fopen(argv[2], "rb");
121
    fseek(f[shift<0], shift < 0 ? -shift : shift, SEEK_SET);
122

    
123
    fseek(f[0],skip_bytes,SEEK_CUR);
124
    fseek(f[1],skip_bytes,SEEK_CUR);
125

    
126
    for(i=0;;){
127
        if( fread(buf[0], SIZE, 1, f[0]) != 1) break;
128
        if( fread(buf[1], SIZE, 1, f[1]) != 1) break;
129

    
130
        for(j=0; j<SIZE; i++,j++){
131
            int64_t a= buf[0][j];
132
            int64_t b= buf[1][j];
133
            if(len==2){
134
                a= (int16_t)(a | (buf[0][++j]<<8));
135
                b= (int16_t)(b | (buf[1][  j]<<8));
136
            }
137
            sse += (a-b) * (a-b);
138
        }
139
    }
140

    
141
    if(!i) i=1;
142
    dev= int_sqrt( ((sse/i)*F*F) + (((sse%i)*F*F) + i/2)/i );
143
    if(sse)
144
        psnr= ((2*log16(max<<16) + log16(i) - log16(sse))*284619LL*F + (1<<31)) / (1LL<<32);
145
    else
146
        psnr= 100*F-1; //floating point free infinity :)
147

    
148
    printf("stddev:%3d.%02d PSNR:%2d.%02d bytes:%d\n",
149
        (int)(dev/F), (int)(dev%F),
150
        (int)(psnr/F), (int)(psnr%F),
151
        i*len);
152
    return 0;
153
}
154

    
155