Statistics
| Branch: | Revision:

ffmpeg / tests / tiny_psnr.c @ 47805d21

History | View | Annotate | Download (3.35 KB)

1 67cbe681 Michael Niedermayer
/*
2
 * Copyright (c) 2003 Michael Niedermayer <michaelni@gmx.at>
3
 *
4
 * This library is free software; you can redistribute it and/or
5
 * modify it under the terms of the GNU Lesser General Public
6
 * License as published by the Free Software Foundation; either
7
 * version 2 of the License, or (at your option) any later version.
8
 *
9
 * This library is distributed in the hope that it will be useful,
10
 * but WITHOUT ANY WARRANTY; without even the implied warranty of
11
 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
12
 * Lesser General Public License for more details.
13
 *
14
 * You should have received a copy of the GNU Lesser General Public
15
 * License along with this library; if not, write to the Free Software
16
 * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA
17
 *
18
 */
19
20
#include <stdio.h>
21 ba96e97f Michael Niedermayer
#include <stdlib.h>
22 67cbe681 Michael Niedermayer
#include <inttypes.h>
23 9e2a16e1 Michael Niedermayer
#include <assert.h>
24 67cbe681 Michael Niedermayer
25
#define F 100
26
#define SIZE 2048
27
28 0c90161f Michael Niedermayer
uint64_t exp16_table[21]={
29 9e2a16e1 Michael Niedermayer
     65537,
30
     65538,
31
     65540,
32
     65544,
33
     65552,
34
     65568,
35
     65600,
36
     65664,
37
     65793,
38
     66050,
39
     66568,
40
     67616,
41
     69763,
42
     74262,
43
     84150,
44
    108051,
45
    178145,
46
    484249,
47
   3578144,
48
 195360063,
49 0c90161f Michael Niedermayer
 582360139072LL,
50 9e2a16e1 Michael Niedermayer
};
51
#if 1
52
// 16.16 fixpoint exp()
53
static unsigned int exp16(unsigned int a){
54
    int i;
55
    int out= 1<<16;
56
57
    for(i=19;i>=0;i--){
58
        if(a&(1<<i))
59
            out= (out*exp16_table[i] + (1<<15))>>16;
60
    }
61
62
    return out;
63
}
64
// 16.16 fixpoint log()
65 0c90161f Michael Niedermayer
static int64_t log16(uint64_t a){
66 9e2a16e1 Michael Niedermayer
    int i;
67
    int out=0;
68 0c90161f Michael Niedermayer
69
    if(a < 1<<16)
70
        return -log16((1LL<<32) / a);
71 9e2a16e1 Michael Niedermayer
    a<<=16;
72
    
73 0c90161f Michael Niedermayer
    for(i=20;i>=0;i--){
74 8176bd1a Michael Niedermayer
        int64_t b= exp16_table[i];
75
        if(a<(b<<16)) continue;
76 9e2a16e1 Michael Niedermayer
        out |= 1<<i;
77 8176bd1a Michael Niedermayer
        a = ((a/b)<<16) + (((a%b)<<16) + b/2)/b;
78 9e2a16e1 Michael Niedermayer
    }
79
    return out;
80
}
81
82
#endif
83 67cbe681 Michael Niedermayer
static uint64_t int_sqrt(uint64_t a)
84
{
85
    uint64_t ret=0;
86
    int s;
87
    uint64_t ret_sq=0;
88
89
    for(s=31; s>=0; s--){
90
        uint64_t b= ret_sq + (1ULL<<(s*2)) + (ret<<s)*2;
91
        if(b<=a){
92
            ret_sq=b;
93
            ret+= 1ULL<<s;
94
        }
95
    }
96
    return ret;
97
}
98
99
int main(int argc,char* argv[]){
100
    int i, j;
101
    uint64_t sse=0;
102
    uint64_t dev;
103
    FILE *f[2];
104
    uint8_t buf[2][SIZE];
105 9e2a16e1 Michael Niedermayer
    uint64_t psnr;
106 ba96e97f Michael Niedermayer
    int len= argc<4 ? 1 : atoi(argv[3]);
107 0c90161f Michael Niedermayer
    int64_t max= (1<<(8*len))-1;
108 ba96e97f Michael Niedermayer
    int shift= argc<5 ? 0 : atoi(argv[4]);
109 67cbe681 Michael Niedermayer
    
110 0c90161f Michael Niedermayer
    if(argc<3){
111 ba96e97f Michael Niedermayer
        printf("tiny_psnr <file1> <file2> [<elem size> [<shift>]]\n");
112 67cbe681 Michael Niedermayer
        return -1;
113
    }
114
    
115 0d6d0cf9 Michael Niedermayer
    f[0]= fopen(argv[1], "rb");
116
    f[1]= fopen(argv[2], "rb");
117 ba96e97f Michael Niedermayer
    fseek(f[shift<0], shift < 0 ? -shift : shift, SEEK_SET);
118 67cbe681 Michael Niedermayer
119
    for(i=0;;){
120
        if( fread(buf[0], SIZE, 1, f[0]) != 1) break;
121
        if( fread(buf[1], SIZE, 1, f[1]) != 1) break;
122
        
123
        for(j=0; j<SIZE; i++,j++){
124 0c90161f Michael Niedermayer
            int64_t a= buf[0][j];
125
            int64_t b= buf[1][j];
126
            if(len==2){
127
                a= (int16_t)(a | (buf[0][++j]<<8));
128
                b= (int16_t)(b | (buf[1][  j]<<8));
129
            }
130 67cbe681 Michael Niedermayer
            sse += (a-b) * (a-b);
131
        }
132
    }
133
    
134 eeaa742c Michael Niedermayer
    if(!i) i=1;
135 0c90161f Michael Niedermayer
    dev= int_sqrt( ((sse/i)*F*F) + (((sse%i)*F*F) + i/2)/i );
136 9e2a16e1 Michael Niedermayer
    if(sse)
137 0c90161f Michael Niedermayer
        psnr= ((2*log16(max<<16) + log16(i) - log16(sse))*284619LL*F + (1<<31)) / (1LL<<32);
138 9e2a16e1 Michael Niedermayer
    else
139
        psnr= 100*F-1; //floating point free infinity :)
140 67cbe681 Michael Niedermayer
    
141 9e2a16e1 Michael Niedermayer
    printf("stddev:%3d.%02d PSNR:%2d.%02d bytes:%d\n", 
142
        (int)(dev/F), (int)(dev%F), 
143
        (int)(psnr/F), (int)(psnr%F),
144
        i);
145 67cbe681 Michael Niedermayer
    return 0;
146
}
147 9e2a16e1 Michael Niedermayer