Libav 0.7.1
|
00001 /* 00002 * Copyright (c) 2003 Michael Niedermayer <michaelni@gmx.at> 00003 * 00004 * This file is part of Libav. 00005 * 00006 * Libav 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 * Libav 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 Libav; 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 // 16.16 fixpoint log() 00056 static int64_t log16(uint64_t a){ 00057 int i; 00058 int out=0; 00059 00060 if(a < 1<<16) 00061 return -log16((1LL<<32) / a); 00062 a<<=16; 00063 00064 for(i=20;i>=0;i--){ 00065 int64_t b= exp16_table[i]; 00066 if(a<(b<<16)) continue; 00067 out |= 1<<i; 00068 a = ((a/b)<<16) + (((a%b)<<16) + b/2)/b; 00069 } 00070 return out; 00071 } 00072 00073 static uint64_t int_sqrt(uint64_t a) 00074 { 00075 uint64_t ret=0; 00076 int s; 00077 uint64_t ret_sq=0; 00078 00079 for(s=31; s>=0; s--){ 00080 uint64_t b= ret_sq + (1ULL<<(s*2)) + (ret<<s)*2; 00081 if(b<=a){ 00082 ret_sq=b; 00083 ret+= 1ULL<<s; 00084 } 00085 } 00086 return ret; 00087 } 00088 00089 int main(int argc,char* argv[]){ 00090 int i, j; 00091 uint64_t sse=0; 00092 uint64_t dev; 00093 FILE *f[2]; 00094 uint8_t buf[2][SIZE]; 00095 uint64_t psnr; 00096 int len= argc<4 ? 1 : atoi(argv[3]); 00097 int64_t max= (1<<(8*len))-1; 00098 int shift= argc<5 ? 0 : atoi(argv[4]); 00099 int skip_bytes = argc<6 ? 0 : atoi(argv[5]); 00100 int size0=0; 00101 int size1=0; 00102 int maxdist = 0; 00103 00104 if(argc<3){ 00105 printf("tiny_psnr <file1> <file2> [<elem size> [<shift> [<skip bytes>]]]\n"); 00106 printf("WAV headers are skipped automatically.\n"); 00107 return 1; 00108 } 00109 00110 f[0]= fopen(argv[1], "rb"); 00111 f[1]= fopen(argv[2], "rb"); 00112 if(!f[0] || !f[1]){ 00113 fprintf(stderr, "Could not open input files.\n"); 00114 return 1; 00115 } 00116 00117 for (i = 0; i < 2; i++) { 00118 uint8_t *p = buf[i]; 00119 if (fread(p, 1, 12, f[i]) != 12) 00120 return 1; 00121 if (!memcmp(p, "RIFF", 4) && 00122 !memcmp(p+8, "WAVE", 4)) { 00123 if (fread(p, 1, 8, f[i]) != 8) 00124 return 1; 00125 while (memcmp(p, "data", 4)) { 00126 int s = p[4] | p[5]<<8 | p[6]<<16 | p[7]<<24; 00127 fseek(f[i], s, SEEK_CUR); 00128 if (fread(p, 1, 8, f[i]) != 8) 00129 return 1; 00130 } 00131 } else { 00132 fseek(f[i], -12, SEEK_CUR); 00133 } 00134 } 00135 00136 fseek(f[shift<0], abs(shift), SEEK_CUR); 00137 00138 fseek(f[0],skip_bytes,SEEK_CUR); 00139 fseek(f[1],skip_bytes,SEEK_CUR); 00140 00141 for(;;){ 00142 int s0= fread(buf[0], 1, SIZE, f[0]); 00143 int s1= fread(buf[1], 1, SIZE, f[1]); 00144 00145 for(j=0; j<FFMIN(s0,s1); j++){ 00146 int64_t a= buf[0][j]; 00147 int64_t b= buf[1][j]; 00148 int dist; 00149 if(len==2){ 00150 a= (int16_t)(a | (buf[0][++j]<<8)); 00151 b= (int16_t)(b | (buf[1][ j]<<8)); 00152 } 00153 sse += (a-b) * (a-b); 00154 dist = abs(a-b); 00155 if (dist > maxdist) maxdist = dist; 00156 } 00157 size0 += s0; 00158 size1 += s1; 00159 if(s0+s1<=0) 00160 break; 00161 } 00162 00163 i= FFMIN(size0,size1)/len; 00164 if(!i) i=1; 00165 dev= int_sqrt( ((sse/i)*F*F) + (((sse%i)*F*F) + i/2)/i ); 00166 if(sse) 00167 psnr= ((2*log16(max<<16) + log16(i) - log16(sse))*284619LL*F + (1LL<<31)) / (1LL<<32); 00168 else 00169 psnr= 1000*F-1; //floating point free infinity :) 00170 00171 printf("stddev:%5d.%02d PSNR:%3d.%02d MAXDIFF:%5d bytes:%9d/%9d\n", 00172 (int)(dev/F), (int)(dev%F), 00173 (int)(psnr/F), (int)(psnr%F), 00174 maxdist, 00175 size0, size1); 00176 return 0; 00177 } 00178 00179