pure sdk for main

This commit is contained in:
divadiow
2025-08-27 09:51:58 +01:00
parent f0d033f1c9
commit 0571416e7c
3283 changed files with 1577720 additions and 1 deletions
+184
View File
@@ -0,0 +1,184 @@
/*Copyright (C) 2008-2009 Timothy B. Terriberry (tterribe@xiph.org)
You can redistribute this library and/or modify it under the terms of the
GNU Lesser General Public License as published by the Free Software
Foundation; either version 2.1 of the License, or (at your option) any later
version.*/
#include "bch15_5.h"
/*A cycle in GF(2**4) generated by alpha=(x**4+x+1).
It is extended an extra 16 entries to avoid some expensive mod operations.*/
static const unsigned char gf16_exp[31]={
1,2,4,8,3,6,12,11,5,10,7,14,15,13,9,1,2,4,8,3,6,12,11,5,10,7,14,15,13,9,1
};
/*The location of each integer 1...16 in the cycle.*/
static const signed char gf16_log[16]={
-1,0,1,4,2,8,5,10,3,14,9,7,6,13,11,12
};
/*Multiplication in GF(2**4) using logarithms.*/
static unsigned gf16_mul(unsigned _a,unsigned _b){
return _a==0||_b==0?0:gf16_exp[gf16_log[_a]+gf16_log[_b]];
}
/*Division in GF(2**4) using logarithms.
The result when dividing by zero is undefined.*/
static unsigned gf16_div(unsigned _a,unsigned _b){
return _a==0?0:gf16_exp[gf16_log[_a]+15-gf16_log[_b]];
}
/*Multiplication in GF(2**4) when the second argument is known to be non-zero
(proven by representing it by its logarithm).*/
static unsigned gf16_hmul(unsigned _a,unsigned _logb){
return _a==0?0:gf16_exp[gf16_log[_a]+_logb];
}
/*The syndrome normally has five values, S_1 ... S_5.
We only calculate and store the odd ones in _s, since S_2=S_1**2 and
S_4=S_2**2.
Returns zero iff all the syndrome values are zero.*/
static int bch15_5_calc_syndrome(unsigned _s[3],unsigned _y){
unsigned p;
int i;
int j;
p=0;
for(i=0;i<15;i++)if(_y&1<<i)p^=gf16_exp[i];
_s[0]=p;
p=0;
for(i=0;i<3;i++)for(j=0;j<5;j++)if(_y&1<<(5*i+j))p^=gf16_exp[j*3];
_s[1]=p;
p=0;
for(i=0;i<5;i++)for(j=0;j<3;j++)if(_y&1<<(3*i+j))p^=gf16_exp[j*5];
_s[2]=p;
return _s[0]!=0||_s[1]!=0||_s[2]!=0;
}
/*Compute the coefficients of the error-locator polynomial.
Returns the number of errors (the degree of the polynomial).*/
static int bch15_5_calc_omega(unsigned _o[3],unsigned _s[3]){
unsigned s02;
unsigned tt;
unsigned dd;
int d;
_o[0]=_s[0];
s02=gf16_mul(_s[0],_s[0]);
dd=_s[1]^gf16_mul(_s[0],s02);
tt=_s[2]^gf16_mul(s02,_s[1]);
_o[1]=dd?gf16_div(tt,dd):0;
_o[2]=dd^gf16_mul(_s[0],_o[1]);
for(d=3;d>0&&!_o[d-1];d--);
return d;
}
/*Find the roots of the error polynomial.
Returns the number of roots found, or a negative value if the polynomial did
not have enough roots, indicating a decoding error.*/
static int bch15_5_calc_epos(unsigned _epos[3],unsigned _s[3]){
unsigned o[3];
int nerrors;
int d;
int i;
d=bch15_5_calc_omega(o,_s);
nerrors=0;
if(d==1)_epos[nerrors++]=gf16_log[o[0]];
else if(d>0){
for(i=0;i<15;i++){
int i2;
i2=gf16_log[gf16_exp[i<<1]];
if(!(gf16_exp[i+i2]^gf16_hmul(o[0],i2)^gf16_hmul(o[1],i)^o[2])){
_epos[nerrors++]=i;
}
}
if(nerrors<d)return -1;
}
return nerrors;
}
int bch15_5_correct(unsigned *_y){
unsigned s[3];
unsigned epos[3];
unsigned y;
int nerrors;
int i;
y=*_y;
if(!bch15_5_calc_syndrome(s,y))return 0;
nerrors=bch15_5_calc_epos(epos,s);
if(nerrors>0){
/*If we had a non-zero syndrome value, we should always find at least one
error location, or we've got a decoding error.*/
for(i=0;i<nerrors;i++)y^=1<<epos[i];
/*If there were too many errors, we may not find enough roots to reduce the
syndrome to zero.
We could recompute it to check, but it's much faster just to check that
we have a valid codeword.*/
if(bch15_5_encode(y>>10)==y){
/*Decoding succeeded.*/
*_y=y;
return nerrors;
}
}
/*Decoding failed due to too many bit errors.*/
return -1;
}
unsigned bch15_5_encode(unsigned _x){
return (-(_x&1)&0x0537)^(-(_x>>1&1)&0x0A6E)^(-(_x>>2&1)&0x11EB)^
(-(_x>>3&1)&0x23D6)^(-(_x>>4&1)&0x429B);
}
#if 0
#include <stdio.h>
static unsigned codes[32];
static int hamming(int _a,int _b){
int d;
int n;
d=_a^_b;
for(n=0;d;n++)d&=d-1;
return n;
}
static int closest(int _y){
int min_i;
int min_d;
int i;
int d;
min_i=0;
min_d=hamming(_y,codes[0]);
for(i=1;i<32;i++){
d=hamming(_y,codes[i]);
if(d<min_d){
min_d=d;
min_i=i;
}
}
return codes[min_i];
}
int main(void){
int i;
/*Print a list of the valid (uncorrupt) codewords.*/
for(i=0;i<32;i++)codes[i]=bch15_5_encode(i);
for(i=0;i<32;i++)printf("0x%04X%s",codes[i],i+1<32?" ":"\n");
/*Try to decode all receivable (possibly corrupt) codewords.*/
for(i=0;i<0x8000;i++){
unsigned y;
unsigned z;
int nerrors;
int j;
y=i;
nerrors=bch15_5_correct(&y);
z=closest(i);
if(nerrors<0){
printf("0x%04X->Failed\n",i);
if(hamming(i,z)<=3)printf("Error: 0x%04X should map to 0x%04X\n",i,z);
}
else{
printf("0x%04X->0x%04X\n",i,y);
if(z!=y)printf("Error: 0x%04X should map to 0x%04X\n",i,z);
}
}
return 0;
}
#endif
+20
View File
@@ -0,0 +1,20 @@
/*Copyright (C) 2008-2009 Timothy B. Terriberry (tterribe@xiph.org)
You can redistribute this library and/or modify it under the terms of the
GNU Lesser General Public License as published by the Free Software
Foundation; either version 2.1 of the License, or (at your option) any later
version.*/
#if !defined(_bch15_5_H)
# define _bch15_5_H (1)
/*Encodes a raw 5-bit value _x into a 15-bit BCH(15,5) code.
This is capable of correcting up to 3 bit errors, and detecting as many as
5 bit errors in some cases.*/
unsigned bch15_5_encode(unsigned _x);
/*Corrects the received code *_y, if possible.
The original data is located in the top five bits.
Returns the number of errors corrected, or a negative value if decoding
failed due to too many bit errors, in which case *_y is left unchanged.*/
int bch15_5_correct(unsigned *_y);
#endif
+740
View File
@@ -0,0 +1,740 @@
/*Copyright (C) 2008-2009 Timothy B. Terriberry (tterribe@xiph.org)
You can redistribute this library and/or modify it under the terms of the
GNU Lesser General Public License as published by the Free Software
Foundation; either version 2.1 of the License, or (at your option) any later
version.*/
#include <stdlib.h>
#include <math.h>
#include <string.h>
#include "util.h"
//#include "type.h"
#include "image.h"
#include "binarize.h"
#if 0
/*Binarization based on~\cite{GPP06}.
@ARTICLE{GPP06,
author="Basilios Gatos and Ioannis E. Pratikakis and Stavros J. Perantonis",
title="Adaptive Degraded Document Image Binarization",
journal="Pattern Recognition",
volume=39,
number=3,
pages="317-327",
month=Mar,
year=2006
}*/
#if 0
/*Applies a 5x5 Wiener filter to the image, in-place, emphasizing differences
where the local variance is small, and de-emphasizing them where it is
large.*/
void qr_wiener_filter(unsigned char *_img,int _width,int _height){
unsigned *m_buf[8];
unsigned *sn2_buf[8];
unsigned char g;
int x;
int y;
if(_width<=0||_height<=0)return;
m_buf[0]=(unsigned *)malloc((_width+4<<3)*sizeof(*m_buf));
sn2_buf[0]=(unsigned *)malloc((_width+4<<3)*sizeof(*sn2_buf));
for(y=1;y<8;y++){
m_buf[y]=m_buf[y-1]+_width+4;
sn2_buf[y]=sn2_buf[y-1]+_width+4;
}
for(y=-4;y<_height;y++){
unsigned *pm;
unsigned *psn2;
int i;
int j;
pm=m_buf[y+2&7];
psn2=sn2_buf[y+2&7];
for(x=-4;x<_width;x++){
unsigned m;
unsigned m2;
m=m2=0;
if(y>=0&&y<_height-4&&x>=0&&x<_width-4)for(i=0;i<5;i++)for(j=0;j<5;j++){
g=_img[(y+i)*_width+x+j];
m+=g;
m2+=g*g;
}
else for(i=0;i<5;i++)for(j=0;j<5;j++){
g=_img[QR_CLAMPI(0,y+i,_height-1)*_width+QR_CLAMPI(0,x+j,_width-1)];
m+=g;
m2+=g*g;
}
pm[x+4]=m;
psn2[x+4]=(m2*25-m*m);
}
pm=m_buf[y&7];
if(y>=0)for(x=0;x<_width;x++){
int sn2;
sn2=sn2_buf[y&7][x+2];
if(sn2){
int vn3;
int m;
/*Gatos et al. give the expression
mu+(s2-v2)*(g-mu)/s2 ,
which we reduce to
mu+(s2-v2)*g/s2-(s2-v2)*mu/s2 ,
g-(v2/s2)*g+(v2/s2)*mu ,
g+(mu-g)*(v2/s2) .
However, s2 is much noisier than v2, and dividing by it often gives
extremely large adjustments, causing speckle near edges.
Therefore we limit the ratio (v2/s2) to lie between 0 and 1.*/
vn3=0;
for(i=-2;i<3;i++){
psn2=sn2_buf[y+i&7];
for(j=0;j<5;j++)vn3+=psn2[x+j];
}
m=m_buf[y&7][x+2];
vn3=vn3+1023>>10;
sn2=25*sn2+1023>>10;
if(vn3<sn2){
int a;
g=_img[y*_width+x];
a=(m-25*g)*vn3;
sn2*=25;
_img[y*_width+x]=QR_CLAMP255(g+QR_DIVROUND(a,sn2));
}
else _img[y*_width+x]=(unsigned char)(((m<<1)+25)/50);
}
}
}
free(sn2_buf[0]);
free(m_buf[0]);
}
#else
/*Applies a 3x3 Wiener filter to the image, in-place, emphasizing differences
where the local variance is small, and de-emphasizing them where it is
large.*/
void qr_wiener_filter(unsigned char *_img,int _width,int _height){
unsigned *m_buf[4];
unsigned *sn2_buf[4];
unsigned char g;
int x;
int y;
if(_width<=0||_height<=0)return;
m_buf[0]=(unsigned *)malloc((_width+2<<2)*sizeof(*m_buf));
sn2_buf[0]=(unsigned *)malloc((_width+2<<2)*sizeof(*sn2_buf));
for(y=1;y<4;y++){
m_buf[y]=m_buf[y-1]+_width+2;
sn2_buf[y]=sn2_buf[y-1]+_width+2;
}
for(y=-2;y<_height;y++){
unsigned *pm;
unsigned *psn2;
int i;
int j;
pm=m_buf[y+1&3];
psn2=sn2_buf[y+1&3];
for(x=-2;x<_width;x++){
unsigned m;
unsigned m2;
m=m2=0;
if(y>=0&&y<_height-2&&x>=0&&x<_width-2)for(i=0;i<3;i++)for(j=0;j<3;j++){
g=_img[(y+i)*_width+x+j];
m+=g;
m2+=g*g;
}
else for(i=0;i<3;i++)for(j=0;j<3;j++){
g=_img[QR_CLAMPI(0,y+i,_height-1)*_width+QR_CLAMPI(0,x+j,_width-1)];
m+=g;
m2+=g*g;
}
pm[x+2]=m;
psn2[x+2]=(m2*9-m*m);
}
pm=m_buf[y&3];
if(y>=0)for(x=0;x<_width;x++){
int sn2;
sn2=sn2_buf[y&3][x+1];
if(sn2){
int m;
int vn3;
/*Gatos et al. give the expression
mu+(s2-v2)*(g-mu)/s2 ,
which we reduce to
mu+(s2-v2)*g/s2-(s2-v2)*mu/s2 ,
g-(v2/s2)*g+(v2/s2)*mu ,
g+(mu-g)*(v2/s2) .
However, s2 is much noisier than v2, and dividing by it often gives
extremely large adjustments, causing speckle near edges.
Therefore we limit the ratio (v2/s2) to lie between 0 and 1.*/
vn3=0;
for(i=-1;i<2;i++){
psn2=sn2_buf[y+i&3];
for(j=0;j<3;j++)vn3+=psn2[x+j];
}
m=m_buf[y&3][x+1];
vn3=vn3+31>>5;
sn2=9*sn2+31>>5;
if(vn3<sn2){
int a;
g=_img[y*_width+x];
a=m-9*g;
sn2*=9;
_img[y*_width+x]=QR_CLAMP255(g+QR_DIVROUND(a,sn2));
}
else _img[y*_width+x]=(unsigned char)(((m<<1)+9)/18);
}
}
}
free(sn2_buf[0]);
free(m_buf[0]);
}
#endif
/*Computes a (conservative) foreground mask using the adaptive binarization
threshold given in~\cite{SP00}, but knocking the threshold parameter down to
k=0.2.
Note on dynamic range: we assume _width*_height<=0x1000000 (24 bits).
Returns the average background value.
@ARTICLE{SP00,
author="Jaakko J. Sauvola and Matti Pietik\"{a}inen",
title="Adaptive Document Image Binarization",
volume=33,
number=2,
pages="225--236",
month=Feb,
year=2000
}*/
static void qr_sauvola_mask(unsigned char *_mask,unsigned *_b,int *_nb,
const unsigned char *_img,int _width,int _height){
unsigned b;
int nb;
b=0;
nb=0;
if(_width>0&&_height>0){
unsigned *col_sums;
unsigned *col2_sums;
int logwindw;
int logwindh;
int windw;
int windh;
int y0offs;
int y1offs;
unsigned g;
unsigned g2;
int x;
int y;
/*We keep the window size fairly large to ensure it doesn't fit completely
inside the center of a finder pattern of a version 1 QR code at full
resolution.*/
for(logwindw=4;logwindw<8&&(1<<logwindw)<(_width+7>>3);logwindw++);
for(logwindh=4;logwindh<8&&(1<<logwindh)<(_height+7>>3);logwindh++);
windw=1<<logwindw;
windh=1<<logwindh;
col_sums=(unsigned *)malloc(_width*sizeof(*col_sums));
col2_sums=(unsigned *)malloc(_width*sizeof(*col2_sums));
/*Initialize sums down each column.*/
for(x=0;x<_width;x++){
g=_img[x];
g2=g*g;
col_sums[x]=(g<<logwindh-1)+g;
col2_sums[x]=(g2<<logwindh-1)+g2;
}
for(y=1;y<(windh>>1);y++){
y1offs=QR_MINI(y,_height-1)*_width;
for(x=0;x<_width;x++){
g=_img[y1offs+x];
col_sums[x]+=g;
col2_sums[x]+=g*g;
}
}
for(y=0;y<_height;y++){
unsigned m;
unsigned m2;
int x0;
int x1;
/*Initialize the sums over the window.*/
m=(col_sums[0]<<logwindw-1)+col_sums[0];
m2=(col2_sums[0]<<logwindw-1)+col2_sums[0];
for(x=1;x<(windw>>1);x++){
x1=QR_MINI(x,_width-1);
m+=col_sums[x1];
m2+=col2_sums[x1];
}
for(x=0;x<_width;x++){
int d;
/*Perform the test against the threshold T = (m/n)*(1+k*(s/R-1)),
where n=windw*windh, s=sqrt((m2-(m*m)/n)/n), and R=128.
We don't actually compute the threshold directly, as that would
require a square root.
Instead we perform the equivalent test:
(m/n)*(m/n)*(m2/n-(m/n)*(m/n))/16 > (((1/k)*g-((1-k)/k)*(m/n))*32)**2
R is split up across each side of the inequality to maximize the
dynamic range available for the right hand side, which requires
31 bits in the worst case.*/
/*(m/n)*(1+(1/5)*(sqrt((m2-m*m/n)/n)/128-1)) > g
m*(1+(1/5)*(sqrt((m2-m*m/n)/n)/128-1)) > g*n
m*sqrt((m2-m*m/n)/n) > 5*g*n-4*m<<7
m*m*(m2*n-m*m) > (5*g*n-4*m<<7)**2*n*n || 5*g*n-4*m < 0 */
g=_img[y*_width+x];
d=(5*g<<logwindw+logwindh)-4*m;
if(d>=0){
unsigned mm;
unsigned mms2;
unsigned d2;
mm=(m>>logwindw)*(m>>logwindh);
mms2=(m2-mm>>logwindw+logwindh)*(mm>>logwindw+logwindh)+15>>4;
d2=d>>logwindw+logwindh-5;
d2*=d2;
if(d2>=mms2){
/*Update the background average.*/
b+=g;
nb++;
_mask[y*_width+x]=0;
}
else _mask[y*_width+x]=0xFF;
}
else _mask[y*_width+x]=0xFF;
/*Update the window sums.*/
if(x+1<_width){
x0=QR_MAXI(0,x-(windw>>1));
x1=QR_MINI(x+(windw>>1),_width-1);
m+=col_sums[x1]-col_sums[x0];
m2+=col2_sums[x1]-col2_sums[x0];
}
}
/*Update the column sums.*/
if(y+1<_height){
y0offs=QR_MAXI(0,y-(windh>>1))*_width;
y1offs=QR_MINI(y+(windh>>1),_height-1)*_width;
for(x=0;x<_width;x++){
g=_img[y0offs+x];
col_sums[x]-=g;
col2_sums[x]-=g*g;
g=_img[y1offs+x];
col_sums[x]+=g;
col2_sums[x]+=g*g;
}
}
}
free(col2_sums);
free(col_sums);
}
*_b=b;
*_nb=nb;
}
/*Interpolates a background image given the source and a conservative
foreground mask.
If the current window contains no foreground pixels, the average background
value over the whole image is used.
Note on dynamic range: we assume _width*_height<=0x8000000 (23 bits).
Returns the average difference between the foreground and the interpolated
background.*/
static void qr_interpolate_background(unsigned char *_dst,
int *_delta,int *_ndelta,const unsigned char *_img,const unsigned char *_mask,
int _width,int _height,unsigned _b,int _nb){
int delta;
int ndelta;
delta=ndelta=0;
if(_width>0&&_height>0){
unsigned *col_sums;
unsigned *ncol_sums;
int logwindw;
int logwindh;
int windw;
int windh;
int y0offs;
int y1offs;
unsigned b;
unsigned g;
int x;
int y;
b=_nb>0?((_b<<1)+_nb)/(_nb<<1):0xFF;
for(logwindw=4;logwindw<8&&(1<<logwindw)<(_width+15>>4);logwindw++);
for(logwindh=4;logwindh<8&&(1<<logwindh)<(_height+15>>4);logwindh++);
windw=1<<logwindw;
windh=1<<logwindh;
col_sums=(unsigned *)malloc(_width*sizeof(*col_sums));
ncol_sums=(unsigned *)malloc(_width*sizeof(*ncol_sums));
/*Initialize sums down each column.*/
for(x=0;x<_width;x++){
if(!_mask[x]){
g=_img[x];
col_sums[x]=(g<<logwindh-1)+g;
ncol_sums[x]=(1<<logwindh-1)+1;
}
else col_sums[x]=ncol_sums[x]=0;
}
for(y=1;y<(windh>>1);y++){
y1offs=QR_MINI(y,_height-1)*_width;
for(x=0;x<_width;x++)if(!_mask[y1offs+x]){
col_sums[x]+=_img[y1offs+x];
ncol_sums[x]++;
}
}
for(y=0;y<_height;y++){
unsigned n;
unsigned m;
int x0;
int x1;
/*Initialize the sums over the window.*/
m=(col_sums[0]<<logwindw-1)+col_sums[0];
n=(ncol_sums[0]<<logwindw-1)+ncol_sums[0];
for(x=1;x<(windw>>1);x++){
x1=QR_MINI(x,_width-1);
m+=col_sums[x1];
n+=ncol_sums[x1];
}
for(x=0;x<_width;x++){
if(!_mask[y*_width+x])g=_img[y*_width+x];
else{
g=n>0?((m<<1)+n)/(n<<1):b;
delta+=(int)g-_img[y*_width+x];
ndelta++;
}
_dst[y*_width+x]=(unsigned char)g;
/*Update the window sums.*/
if(x+1<_width){
x0=QR_MAXI(0,x-(windw>>1));
x1=QR_MINI(x+(windw>>1),_width-1);
m+=col_sums[x1]-col_sums[x0];
n+=ncol_sums[x1]-ncol_sums[x0];
}
}
/*Update the column sums.*/
if(y+1<_height){
y0offs=QR_MAXI(0,y-(windh>>1))*_width;
y1offs=QR_MINI(y+(windh>>1),_height-1)*_width;
for(x=0;x<_width;x++){
if(!_mask[y0offs+x]){
col_sums[x]-=_img[y0offs+x];
ncol_sums[x]--;
}
if(!_mask[y1offs+x]){
col_sums[x]+=_img[y1offs+x];
ncol_sums[x]++;
}
}
}
}
free(ncol_sums);
free(col_sums);
}
*_delta=delta;
*_ndelta=ndelta;
}
/*Parameters of the logistic sigmoid function that defines the threshold based
on the background intensity.
They should all be between 0 and 1.*/
#define QR_GATOS_Q (0.7)
#define QR_GATOS_P1 (0.5)
#define QR_GATOS_P2 (0.8)
/*Compute the final binarization mask according to Gatos et al.'s
method~\cite{GPP06}.*/
static void qr_gatos_mask(unsigned char *_mask,const unsigned char *_img,
const unsigned char *_background,int _width,int _height,
unsigned _b,int _nb,int _delta,int _ndelta){
unsigned thresh[256];
unsigned g;
double delta;
double b;
int x;
int y;
/*Construct a lookup table for the thresholds.
This bit uses floating point, but doesn't need to do much calculation, so
emulation should be fine.*/
b=_nb>0?(_b+0.5)/_nb:0xFF;
delta=_ndelta>0?(_delta+0.5)/_ndelta:0xFF;
for(g=0;g<256;g++){
double d;
d=QR_GATOS_Q*delta*(QR_GATOS_P2+(1-QR_GATOS_P2)/
(1+exp(2*(1+QR_GATOS_P1)/(1-QR_GATOS_P1)-4*g/(b*(1-QR_GATOS_P1)))));
if(d<1)d=1;
else if(d>0xFF)d=0xFF;
thresh[g]=(unsigned)floor(d);
}
/*Apply the adaptive threshold.*/
for(y=0;y<_height;y++)for(x=0;x<_width;x++){
g=_background[y*_width+x];
/*_background[y*_width+x]=thresh[g];*/
_mask[y*_width+x]=(unsigned char)(-(g-_img[y*_width+x]>thresh[g])&0xFF);
}
/*{
FILE *fout;
fout=fopen("thresh.png","wb");
image_write_png(_background,_width,_height,fout);
fclose(fout);
}*/
}
/*Binarizes a grayscale image.*/
void qr_binarize(unsigned char *_img,int _width,int _height){
unsigned char *mask;
unsigned char *background;
unsigned b;
int nb;
int delta;
int ndelta;
/*qr_wiener_filter(_img,_width,_height);
{
FILE *fout;
fout=fopen("wiener.png","wb");
image_write_png(_img,_width,_height,fout);
fclose(fout);
}*/
mask=(unsigned char *)malloc(_width*_height*sizeof(*mask));
qr_sauvola_mask(mask,&b,&nb,_img,_width,_height);
/*{
FILE *fout;
fout=fopen("foreground.png","wb");
image_write_png(mask,_width,_height,fout);
fclose(fout);
}*/
background=(unsigned char *)malloc(_width*_height*sizeof(*mask));
qr_interpolate_background(background,&delta,&ndelta,
_img,mask,_width,_height,b,nb);
/*{
FILE *fout;
fout=fopen("background.png","wb");
image_write_png(background,_width,_height,fout);
fclose(fout);
}*/
qr_gatos_mask(_img,_img,background,_width,_height,b,nb,delta,ndelta);
free(background);
free(mask);
}
#else
/*The above algorithms are computationally expensive, and do not work as well
as the simple algorithm below.
Sauvola by itself does an excellent job of classifying regions outside the
QR code as background, which greatly reduces the chance of false alarms.
However, it also tends to over-shrink isolated black dots inside the code,
making them easy to miss with even slight mis-alignment.
Since the Gatos method uses Sauvola as input to its background interpolation
method, it cannot possibly mark any pixels as foreground which Sauvola
classified as background, and thus suffers from the same problem.
The following simple adaptive threshold method does not have this problem,
though it produces essentially random noise outside the QR code region.
QR codes are structured well enough that this does not seem to lead to any
actual false alarms in practice, and it allows many more codes to be
detected and decoded successfully than the Sauvola or Gatos binarization
methods.*/
/*A simplified adaptive thresholder.
This compares the current pixel value to the mean value of a (large) window
surrounding it.*/
#if 0
unsigned char psram_bin[640*480]__attribute__ ((section(".psram.src")));
unsigned char *qr_binarize(const unsigned char *_img,int _width,int _height){
unsigned char *mask = NULL;
if(_width>0&&_height>0){
unsigned *col_sums;
int logwindw;
int logwindh;
int windw;
int windh;
int y0offs;
int y1offs;
unsigned g;
int x;
int y;
//mask=(unsigned char *)malloc(_width*_height*sizeof(*mask));
mask = psram_bin;
/*We keep the window size fairly large to ensure it doesn't fit completely
inside the center of a finder pattern of a version 1 QR code at full
resolution.*/
for(logwindw=4;logwindw<8&&(1<<logwindw)<(_width+7>>3);logwindw++);
for(logwindh=4;logwindh<8&&(1<<logwindh)<(_height+7>>3);logwindh++);
windw=1<<logwindw;
windh=1<<logwindh;
col_sums=(unsigned *)malloc(_width*sizeof(*col_sums));
/*Initialize sums down each column.*/
for(x=0;x<_width;x++){
g=_img[x];
col_sums[x]=(g<<logwindh-1)+g;
}
for(y=1;y<(windh>>1);y++){
y1offs=QR_MINI(y,_height-1)*_width;
for(x=0;x<_width;x++){
g=_img[y1offs+x];
col_sums[x]+=g;
}
}
for(y=0;y<_height;y++){
unsigned m;
int x0;
int x1;
/*Initialize the sum over the window.*/
m=(col_sums[0]<<logwindw-1)+col_sums[0];
for(x=1;x<(windw>>1);x++){
x1=QR_MINI(x,_width-1);
m+=col_sums[x1];
}
for(x=0;x<_width;x++){
/*Perform the test against the threshold T = (m/n)-D,
where n=windw*windh and D=3.*/
g=_img[y*_width+x];
mask[y*_width+x]=-(g+3<<logwindw+logwindh<m)&0xFF;
//printf("g:%x m:%x logwindh:%x,logwindw:%x,_width:%x,x:%x,y:%x,mask[y*_width+x]:%x\r\n",g,m,logwindh,logwindw,_width,x,y,mask[y*_width+x]);
/*Update the window sum.*/
if(x+1<_width){
x0=QR_MAXI(0,x-(windw>>1));
x1=QR_MINI(x+(windw>>1),_width-1);
m+=col_sums[x1]-col_sums[x0];
}
}
/*Update the column sums.*/
if(y+1<_height){
y0offs=QR_MAXI(0,y-(windh>>1))*_width;
y1offs=QR_MINI(y+(windh>>1),_height-1)*_width;
for(x=0;x<_width;x++){
col_sums[x]-=_img[y0offs+x];
col_sums[x]+=_img[y1offs+x];
}
}
}
free(col_sums);
}
#if defined(QR_DEBUG)
{
FILE *fout;
fout=fopen("binary.png","wb");
image_write_png(_img,_width,_height,fout);
fclose(fout);
}
#endif
return(mask);
}
#endif
//与qr_binarize一致,只是mask的空间从外部申请然后匹配
unsigned char *qr_binarize2(unsigned char *mask,const unsigned char *_img,int _width,int _height){
if(_width>0&&_height>0){
unsigned *col_sums;
int logwindw;
int logwindh;
int windw;
int windh;
int y0offs;
int y1offs;
unsigned g;
int x;
int y;
/*We keep the window size fairly large to ensure it doesn't fit completely
inside the center of a finder pattern of a version 1 QR code at full
resolution.*/
for(logwindw=4;logwindw<8&&(1<<logwindw)<((_width+7)>>3);logwindw++);
for(logwindh=4;logwindh<8&&(1<<logwindh)<((_height+7)>>3);logwindh++);
windw=1<<logwindw;
windh=1<<logwindh;
col_sums=(unsigned *)malloc(_width*sizeof(*col_sums));
/*Initialize sums down each column.*/
for(x=0;x<_width;x++){
g=_img[x];
col_sums[x]=(g<<(logwindh-1))+g;
}
for(y=1;y<(windh>>1);y++){
y1offs=QR_MINI(y,_height-1)*_width;
for(x=0;x<_width;x++){
g=_img[y1offs+x];
col_sums[x]+=g;
}
}
for(y=0;y<_height;y++){
unsigned m;
int x0;
int x1;
/*Initialize the sum over the window.*/
m=(col_sums[0]<<(logwindw-1))+col_sums[0];
for(x=1;x<(windw>>1);x++){
x1=QR_MINI(x,_width-1);
m+=col_sums[x1];
}
for(x=0;x<_width;x++){
/*Perform the test against the threshold T = (m/n)-D,
where n=windw*windh and D=3.*/
g=_img[y*_width+x];
mask[y*_width+x]=-((g+3)<<(logwindw+logwindh)<m)&0xFF;
//printf("g:%x m:%x logwindh:%x,logwindw:%x,_width:%x,x:%x,y:%x,mask[y*_width+x]:%x\r\n",g,m,logwindh,logwindw,_width,x,y,mask[y*_width+x]);
/*Update the window sum.*/
if(x+1<_width){
x0=QR_MAXI(0,x-(windw>>1));
x1=QR_MINI(x+(windw>>1),_width-1);
m+=col_sums[x1]-col_sums[x0];
}
}
/*Update the column sums.*/
if(y+1<_height){
y0offs=QR_MAXI(0,y-(windh>>1))*_width;
y1offs=QR_MINI(y+(windh>>1),_height-1)*_width;
for(x=0;x<_width;x++){
col_sums[x]-=_img[y0offs+x];
col_sums[x]+=_img[y1offs+x];
}
}
}
free(col_sums);
}
#if defined(QR_DEBUG)
{
FILE *fout;
fout=fopen("binary.png","wb");
image_write_png(_img,_width,_height,fout);
fclose(fout);
}
#endif
return(mask);
}
#endif
#if defined(TEST_BINARIZE)
#include <stdio.h>
#include "image.c"
int main(int _argc,char **_argv){
unsigned char *img;
int width;
int height;
int x;
int y;
if(_argc<2){
fprintf(stderr,"usage: %s <image>.png\n",_argv[0]);
return EXIT_FAILURE;
}
/*width=1182;
height=1181;
img=(unsigned char *)malloc(width*height*sizeof(*img));
for(y=0;y<height;y++)for(x=0;x<width;x++){
img[y*width+x]=(unsigned char)(-((x&1)^(y&1))&0xFF);
}*/
{
FILE *fin;
fin=fopen(_argv[1],"rb");
image_read_png(&img,&width,&height,fin);
fclose(fin);
}
qr_binarize(img,width,height);
/*{
FILE *fout;
fout=fopen("binary.png","wb");
image_write_png(img,width,height,fout);
fclose(fout);
}*/
free(img);
return EXIT_SUCCESS;
}
#endif
+17
View File
@@ -0,0 +1,17 @@
/*Copyright (C) 2008-2009 Timothy B. Terriberry (tterribe@xiph.org)
You can redistribute this library and/or modify it under the terms of the
GNU Lesser General Public License as published by the Free Software
Foundation; either version 2.1 of the License, or (at your option) any later
version.*/
#if !defined(_qrcode_binarize_H)
# define _qrcode_binarize_H (1)
void qr_image_cross_masking_median_filter(unsigned char *_img,
int _width,int _height);
void qr_wiener_filter(unsigned char *_img,int _width,int _height);
/*Binarizes a grayscale image.*/
unsigned char *qr_binarize(const unsigned char *_img,int _width,int _height);
#endif
+139
View File
@@ -0,0 +1,139 @@
/*Written by Timothy B. Terriberry (tterribe@xiph.org) 1999-2009 public domain.
Based on the public domain implementation by Robert J. Jenkins Jr.*/
#include <float.h>
#include <math.h>
#include <string.h>
#include "isaac.h"
#define ISAAC_MASK (0xFFFFFFFFU)
static void isaac_update(isaac_ctx *_ctx){
unsigned *m;
unsigned *r;
unsigned a;
unsigned b;
unsigned x;
unsigned y;
int i;
m=_ctx->m;
r=_ctx->r;
a=_ctx->a;
b=(_ctx->b+(++_ctx->c))&ISAAC_MASK;
for(i=0;i<ISAAC_SZ/2;i++){
x=m[i];
a=((a^a<<13)+m[i+ISAAC_SZ/2])&ISAAC_MASK;
m[i]=y=(m[(x&(ISAAC_SZ-1)<<2)>>2]+a+b)&ISAAC_MASK;
r[i]=b=(m[y>>(ISAAC_SZ_LOG+2)&(ISAAC_SZ-1)]+x)&ISAAC_MASK;
x=m[++i];
a=((a^a>>6)+m[i+ISAAC_SZ/2])&ISAAC_MASK;
m[i]=y=(m[(x&(ISAAC_SZ-1)<<2)>>2]+a+b)&ISAAC_MASK;
r[i]=b=(m[y>>(ISAAC_SZ_LOG+2)&(ISAAC_SZ-1)]+x)&ISAAC_MASK;
x=m[++i];
a=((a^a<<2)+m[i+ISAAC_SZ/2])&ISAAC_MASK;
m[i]=y=(m[(x&(ISAAC_SZ-1)<<2)>>2]+a+b)&ISAAC_MASK;
r[i]=b=(m[y>>(ISAAC_SZ_LOG+2)&(ISAAC_SZ-1)]+x)&ISAAC_MASK;
x=m[++i];
a=((a^a>>16)+m[i+ISAAC_SZ/2])&ISAAC_MASK;
m[i]=y=(m[(x&(ISAAC_SZ-1)<<2)>>2]+a+b)&ISAAC_MASK;
r[i]=b=(m[y>>(ISAAC_SZ_LOG+2)&(ISAAC_SZ-1)]+x)&ISAAC_MASK;
}
for(i=ISAAC_SZ/2;i<ISAAC_SZ;i++){
x=m[i];
a=((a^a<<13)+m[i-ISAAC_SZ/2])&ISAAC_MASK;
m[i]=y=(m[(x&(ISAAC_SZ-1)<<2)>>2]+a+b)&ISAAC_MASK;
r[i]=b=(m[y>>(ISAAC_SZ_LOG+2)&(ISAAC_SZ-1)]+x)&ISAAC_MASK;
x=m[++i];
a=((a^a>>6)+m[i-ISAAC_SZ/2])&ISAAC_MASK;
m[i]=y=(m[(x&(ISAAC_SZ-1)<<2)>>2]+a+b)&ISAAC_MASK;
r[i]=b=(m[y>>(ISAAC_SZ_LOG+2)&(ISAAC_SZ-1)]+x)&ISAAC_MASK;
x=m[++i];
a=((a^a<<2)+m[i-ISAAC_SZ/2])&ISAAC_MASK;
m[i]=y=(m[(x&(ISAAC_SZ-1)<<2)>>2]+a+b)&ISAAC_MASK;
r[i]=b=(m[y>>(ISAAC_SZ_LOG+2)&(ISAAC_SZ-1)]+x)&ISAAC_MASK;
x=m[++i];
a=((a^a>>16)+m[i-ISAAC_SZ/2])&ISAAC_MASK;
m[i]=y=(m[(x&(ISAAC_SZ-1)<<2)>>2]+a+b)&ISAAC_MASK;
r[i]=b=(m[y>>(ISAAC_SZ_LOG+2)&(ISAAC_SZ-1)]+x)&ISAAC_MASK;
}
_ctx->b=b;
_ctx->a=a;
_ctx->n=ISAAC_SZ;
}
static void isaac_mix(unsigned _x[8]){
static const unsigned char SHIFT[8]={11,2,8,16,10,4,8,9};
int i;
for(i=0;i<8;i++){
_x[i]^=_x[(i+1)&7]<<SHIFT[i];
_x[(i+3)&7]+=_x[i];
_x[(i+1)&7]+=_x[(i+2)&7];
i++;
_x[i]^=_x[(i+1)&7]>>SHIFT[i];
_x[(i+3)&7]+=_x[i];
_x[(i+1)&7]+=_x[(i+2)&7];
}
}
void isaac_init(isaac_ctx *_ctx,const void *_seed,int _nseed){
const unsigned char *seed;
unsigned *m;
unsigned *r;
unsigned x[8];
int i;
int j;
_ctx->a=_ctx->b=_ctx->c=0;
m=_ctx->m;
r=_ctx->r;
x[0]=x[1]=x[2]=x[3]=x[4]=x[5]=x[6]=x[7]=0x9E3779B9;
for(i=0;i<4;i++)isaac_mix(x);
if(_nseed>ISAAC_SEED_SZ_MAX)_nseed=ISAAC_SEED_SZ_MAX;
seed=(const unsigned char *)_seed;
for(i=0;i<_nseed>>2;i++){
r[i]=seed[i<<2|3]<<24|seed[i<<2|2]<<16|seed[i<<2|1]<<8|seed[i<<2];
}
if(_nseed&3){
r[i]=seed[i<<2];
for(j=1;j<(_nseed&3);j++)r[i]+=seed[i<<2|j]<<(j<<3);
i++;
}
memset(r+i,0,(ISAAC_SZ-i)*sizeof(*r));
for(i=0;i<ISAAC_SZ;i+=8){
for(j=0;j<8;j++)x[j]+=r[i+j];
isaac_mix(x);
memcpy(m+i,x,sizeof(x));
}
for(i=0;i<ISAAC_SZ;i+=8){
for(j=0;j<8;j++)x[j]+=m[i+j];
isaac_mix(x);
memcpy(m+i,x,sizeof(x));
}
isaac_update(_ctx);
}
unsigned isaac_next_uint32(isaac_ctx *_ctx){
if(!_ctx->n)isaac_update(_ctx);
return _ctx->r[--_ctx->n];
}
/*Returns a uniform random integer less than the given maximum value.
_n: The upper bound on the range of numbers returned (not inclusive).
This must be strictly less than 2**32.
Return: An integer uniformly distributed between 0 (inclusive) and _n
(exclusive).*/
unsigned isaac_next_uint(isaac_ctx *_ctx,unsigned _n){
unsigned r;
unsigned v;
unsigned d;
do{
r=isaac_next_uint32(_ctx);
v=r%_n;
d=r-v;
}
while(((d+_n-1)&ISAAC_MASK)<d);
return v;
}
+41
View File
@@ -0,0 +1,41 @@
/*Written by Timothy B. Terriberry (tterribe@xiph.org) 1999-2009 public domain.
Based on the public domain implementation by Robert J. Jenkins Jr.*/
#if !defined(_isaac_H)
# define _isaac_H (1)
typedef struct isaac_ctx isaac_ctx;
#define ISAAC_SZ_LOG (8)
#define ISAAC_SZ (1<<ISAAC_SZ_LOG)
#define ISAAC_SEED_SZ_MAX (ISAAC_SZ<<2)
/*ISAAC is the most advanced of a series of Pseudo-Random Number Generators
designed by Robert J. Jenkins Jr. in 1996.
http://www.burtleburtle.net/bob/rand/isaac.html
To quote:
No efficient method is known for deducing their internal states.
ISAAC requires an amortized 18.75 instructions to produce a 32-bit value.
There are no cycles in ISAAC shorter than 2**40 values.
The expected cycle length is 2**8295 values.*/
struct isaac_ctx{
unsigned n;
unsigned r[ISAAC_SZ];
unsigned m[ISAAC_SZ];
unsigned a;
unsigned b;
unsigned c;
};
void isaac_init(isaac_ctx *_ctx,const void *_seed,int _nseed);
unsigned isaac_next_uint32(isaac_ctx *_ctx);
unsigned isaac_next_uint(isaac_ctx *_ctx,unsigned _n);
#endif
File diff suppressed because it is too large Load Diff
+168
View File
@@ -0,0 +1,168 @@
/*Copyright (C) 2008-2009 Timothy B. Terriberry (tterribe@xiph.org)
You can redistribute this library and/or modify it under the terms of the
GNU Lesser General Public License as published by the Free Software
Foundation; either version 2.1 of the License, or (at your option) any later
version.*/
#if !defined(_qrdec_H)
# define _qrdec_H (1)
#include <zbar.h>
typedef struct qr_code_data_entry qr_code_data_entry;
typedef struct qr_code_data qr_code_data;
typedef struct qr_code_data_list qr_code_data_list;
typedef enum qr_mode{
/*Numeric digits ('0'...'9').*/
QR_MODE_NUM=1,
/*Alphanumeric characters ('0'...'9', 'A'...'Z', plus the punctuation
' ', '$', '%', '*', '+', '-', '.', '/', ':').*/
QR_MODE_ALNUM,
/*Structured-append header.*/
QR_MODE_STRUCT,
/*Raw 8-bit bytes.*/
QR_MODE_BYTE,
/*FNC1 marker (for more info, see http://www.mecsw.com/specs/uccean128.html).
In the "first position" data is formatted in accordance with GS1 General
Specifications.*/
QR_MODE_FNC1_1ST,
/*Mode 6 reserved?*/
/*Extended Channel Interpretation code.*/
QR_MODE_ECI=7,
/*SJIS kanji characters.*/
QR_MODE_KANJI,
/*FNC1 marker (for more info, see http://www.mecsw.com/specs/uccean128.html).
In the "second position" data is formatted in accordance with an industry
application as specified by AIM Inc.*/
QR_MODE_FNC1_2ND
}qr_mode;
/*Check if a mode has a data buffer associated with it.
Currently this is only modes with exactly one bit set.*/
#define QR_MODE_HAS_DATA(_mode) (!((_mode)&((_mode)-1)))
/*ECI may be used to signal a character encoding for the data.*/
typedef enum qr_eci_encoding{
/*GLI0 is like CP437, but the encoding is reset at the beginning of each
structured append symbol.*/
QR_ECI_GLI0,
/*GLI1 is like ISO8859_1, but the encoding is reset at the beginning of each
structured append symbol.*/
QR_ECI_GLI1,
/*The remaining encodings do not reset at the start of the next structured
append symbol.*/
QR_ECI_CP437,
/*Western European.*/
QR_ECI_ISO8859_1,
/*Central European.*/
QR_ECI_ISO8859_2,
/*South European.*/
QR_ECI_ISO8859_3,
/*North European.*/
QR_ECI_ISO8859_4,
/*Cyrillic.*/
QR_ECI_ISO8859_5,
/*Arabic.*/
QR_ECI_ISO8859_6,
/*Greek.*/
QR_ECI_ISO8859_7,
/*Hebrew.*/
QR_ECI_ISO8859_8,
/*Turkish.*/
QR_ECI_ISO8859_9,
/*Nordic.*/
QR_ECI_ISO8859_10,
/*Thai.*/
QR_ECI_ISO8859_11,
/*There is no ISO/IEC 8859-12.*/
/*Baltic rim.*/
QR_ECI_ISO8859_13=QR_ECI_ISO8859_11+2,
/*Celtic.*/
QR_ECI_ISO8859_14,
/*Western European with euro.*/
QR_ECI_ISO8859_15,
/*South-Eastern European (with euro).*/
QR_ECI_ISO8859_16,
/*ECI 000019 is reserved?*/
/*Shift-JIS.*/
QR_ECI_SJIS=20
}qr_eci_encoding;
/*A single unit of parsed QR code data.*/
struct qr_code_data_entry{
/*The mode of this data block.*/
qr_mode mode;
union{
/*Data buffer for modes that have one.*/
struct{
unsigned char *buf;
int len;
}data;
/*Decoded "Extended Channel Interpretation" data.*/
unsigned eci;
/*Structured-append header data.*/
struct{
unsigned char sa_index;
unsigned char sa_size;
unsigned char sa_parity;
}sa;
}payload;
};
/*Low-level QR code data.*/
struct qr_code_data{
/*The decoded data entries.*/
qr_code_data_entry *entries;
int nentries;
/*The code version (1...40).*/
unsigned char version;
/*The ECC level (0...3, corresponding to 'L', 'M', 'Q', and 'H').*/
unsigned char ecc_level;
/*Structured-append information.*/
/*The index of this code in the structured-append group.
If sa_size is zero, this is undefined.*/
unsigned char sa_index;
/*The size of the structured-append group, or 0 if there was no S-A header.*/
unsigned char sa_size;
/*The parity of the entire structured-append group.
If sa_size is zero, this is undefined.*/
unsigned char sa_parity;
/*The parity of this code.
If sa_size is zero, this is undefined.*/
unsigned char self_parity;
/*An approximate bounding box for the code.
Points appear in the order up-left, up-right, down-left, down-right,
relative to the orientation of the QR code.*/
qr_point bbox[4];
};
struct qr_code_data_list{
qr_code_data *qrdata;
int nqrdata;
int cqrdata;
};
/*Extract symbol data from a list of QR codes and attach to the image.
All text is converted to UTF-8.
Any structured-append group that does not have all of its members is decoded
as ZBAR_PARTIAL with ZBAR_PARTIAL components for the discontinuities.
Note that isolated members of a structured-append group may be decoded with
the wrong character set, since the correct setting cannot be propagated
between codes.
Return: The number of symbols which were successfully extracted from the
codes; this will be at most the number of codes.*/
int qr_code_data_list_extract_text(const qr_code_data_list *_qrlist,
zbar_image_scanner_t *iscn,
zbar_image_t *img);
/*TODO: Parse DoCoMo standard barcode data formats.
See http://www.nttdocomo.co.jp/english/service/imode/make/content/barcode/function/application/
for details.*/
#endif
+398
View File
@@ -0,0 +1,398 @@
/*Copyright (C) 2008-2009 Timothy B. Terriberry (tterribe@xiph.org)
You can redistribute this library and/or modify it under the terms of the
GNU Lesser General Public License as published by the Free Software
Foundation; either version 2.1 of the License, or (at your option) any later
version.*/
#include <stdio.h>
#include <stdlib.h>
#include <string.h>
//#include <iconv.h>
#include "qrcode.h"
#include "qrdec.h"
#include "util.h"
//#include "type.h"
#include "image.h"
#include "error.h"
#include "img_scanner.h"
static int text_is_ascii(const unsigned char *_text,int _len){
int i;
for(i=0;i<_len;i++)if(_text[i]>=0x80)return 0;
return 1;
}
static int text_is_latin1(const unsigned char *_text,int _len){
int i;
for(i=0;i<_len;i++){
/*The following line fails to compile correctly with gcc 3.4.4 on ARM with
any optimizations enabled.*/
if(_text[i]>=0x80&&_text[i]<0xA0)return 0;
}
return 1;
}
static void enc_list_mtf(iconv_t _enc_list[3],iconv_t _enc){
int i;
for(i=0;i<3;i++)if(_enc_list[i]==_enc){
int j;
for(j=i;j-->0;)_enc_list[j+1]=_enc_list[j];
_enc_list[0]=_enc;
break;
}
}
extern size_t zbar_iconv(unsigned char* _cd, char **inbuf, size_t *inbytesleft, char **outbuf, size_t *outbytesleft);
extern iconv_t zbar_iconv_open(const char *tocode, const char *fromcode);
extern int zbar_iconv_close(iconv_t _cd);
int qr_code_data_list_extract_text(const qr_code_data_list *_qrlist,
zbar_image_scanner_t *iscn,
zbar_image_t *img)
{
iconv_t sjis_cd;
iconv_t utf8_cd;
iconv_t latin1_cd;
const qr_code_data *qrdata;
int nqrdata;
unsigned char *mark;
char **text;
int ntext;
int i;
qrdata=_qrlist->qrdata;
nqrdata=_qrlist->nqrdata;
text=(char **)malloc(nqrdata*sizeof(*text));
mark=(unsigned char *)calloc(nqrdata,sizeof(*mark));
ntext=0;
/*This is the encoding the standard says is the default.*/
latin1_cd=zbar_iconv_open("UTF-8","ISO8859-1");
/*But this one is often used, as well.*/
sjis_cd=zbar_iconv_open("UTF-8","SJIS");
/*This is a trivial conversion just to check validity without extra code.*/
utf8_cd=zbar_iconv_open("UTF-8","UTF-8");
for(i=0;i<nqrdata;i++)if(!mark[i]){
const qr_code_data *qrdataj;
const qr_code_data_entry *entry;
iconv_t enc_list[3];
iconv_t eci_cd;
int sa[16];
int sa_size;
char *sa_text;
size_t sa_ntext;
size_t sa_ctext;
int fnc1;
int eci;
int err;
int j;
int k;
/*Step 0: Collect the other QR codes belonging to this S-A group.*/
if(qrdata[i].sa_size){
unsigned sa_parity;
sa_size=qrdata[i].sa_size;
sa_parity=qrdata[i].sa_parity;
for(j=0;j<sa_size;j++)sa[j]=-1;
for(j=i;j<nqrdata;j++)if(!mark[j]){
/*TODO: We could also match version, ECC level, etc. if size and
parity alone are too ambiguous.*/
if(qrdata[j].sa_size==sa_size&&qrdata[j].sa_parity==sa_parity&&
sa[qrdata[j].sa_index]<0){
sa[qrdata[j].sa_index]=j;
mark[j]=1;
}
}
/*TODO: If the S-A group is complete, check the parity.*/
}
else{
sa[0]=i;
sa_size=1;
}
sa_ctext=0;
fnc1=0;
/*Step 1: Detect FNC1 markers and estimate the required buffer size.*/
for(j=0;j<sa_size;j++)if(sa[j]>=0){
qrdataj=qrdata+sa[j];
for(k=0;k<qrdataj->nentries;k++){
int shift;
entry=qrdataj->entries+k;
shift=0;
switch(entry->mode){
/*FNC1 applies to the entire code and ignores subsequent markers.*/
case QR_MODE_FNC1_1ST:
case QR_MODE_FNC1_2ND:fnc1=1;break;
/*2 SJIS bytes will be at most 4 UTF-8 bytes.*/
case QR_MODE_KANJI:shift++;
/*We assume at most 4 UTF-8 bytes per input byte.
I believe this is true for all the encodings we actually use.*/
case QR_MODE_BYTE:shift++;
default:{
/*The remaining two modes are already valid UTF-8.*/
if(QR_MODE_HAS_DATA(entry->mode)){
sa_ctext+=entry->payload.data.len<<shift;
}
}break;
}
}
}
/*Step 2: Convert the entries.*/
sa_text=(char *)malloc((sa_ctext+1)*sizeof(*sa_text));
sa_ntext=0;
eci=-1;
enc_list[0]=sjis_cd;
enc_list[1]=latin1_cd;
enc_list[2]=utf8_cd;
eci_cd=(iconv_t)-1;
err=0;
zbar_symbol_t *syms = NULL, **sym = &syms;
for(j = 0; j < sa_size && !err; j++, sym = &(*sym)->next) {
*sym = _zbar_image_scanner_alloc_sym(iscn, ZBAR_QRCODE, 0);
(*sym)->datalen = sa_ntext;
if(sa[j]<0){
/* generic placeholder for unfinished results */
(*sym)->type = ZBAR_PARTIAL;
/*Skip all contiguous missing segments.*/
for(j++;j<sa_size&&sa[j]<0;j++);
/*If there aren't any more, stop.*/
if(j>=sa_size)break;
/* mark break in data */
sa_text[sa_ntext++]='\0';
(*sym)->datalen = sa_ntext;
/* advance to next symbol */
sym = &(*sym)->next;
*sym = _zbar_image_scanner_alloc_sym(iscn, ZBAR_QRCODE, 0);
}
qrdataj=qrdata+sa[j];
/* expose bounding box */
sym_add_point(*sym, qrdataj->bbox[0][0], qrdataj->bbox[0][1]);
sym_add_point(*sym, qrdataj->bbox[2][0], qrdataj->bbox[2][1]);
sym_add_point(*sym, qrdataj->bbox[3][0], qrdataj->bbox[3][1]);
sym_add_point(*sym, qrdataj->bbox[1][0], qrdataj->bbox[1][1]);
for(k=0;k<qrdataj->nentries&&!err;k++){
size_t inleft;
size_t outleft;
char *in;
char *out;
entry=qrdataj->entries+k;
switch(entry->mode){
case QR_MODE_NUM:{
if(sa_ctext-sa_ntext>=(size_t)entry->payload.data.len){
memcpy(sa_text+sa_ntext,entry->payload.data.buf,
entry->payload.data.len*sizeof(*sa_text));
sa_ntext+=entry->payload.data.len;
}
else err=1;
}break;
case QR_MODE_ALNUM:{
char *p;
in=(char *)entry->payload.data.buf;
inleft=entry->payload.data.len;
/*FNC1 uses '%' as an escape character.*/
if(fnc1)for(;;){
size_t plen;
char c;
p=memchr(in,'%',inleft*sizeof(*in));
if(p==NULL)break;
plen=p-in;
if(sa_ctext-sa_ntext<plen+1)break;
memcpy(sa_text+sa_ntext,in,plen*sizeof(*in));
sa_ntext+=plen;
/*Two '%'s is a literal '%'*/
if(plen+1<inleft&&p[1]=='%'){
c='%';
plen++;
p++;
}
/*One '%' is the ASCII group separator.*/
else c=0x1D;
sa_text[sa_ntext++]=c;
inleft-=plen+1;
in=p+1;
}
else p=NULL;
if(p!=NULL||sa_ctext-sa_ntext<inleft)err=1;
else{
memcpy(sa_text+sa_ntext,in,inleft*sizeof(*sa_text));
sa_ntext+=inleft;
}
}break;
/*TODO: This will not handle a multi-byte sequence split between
multiple data blocks.
Does such a thing occur?
Is it allowed?
It requires copying buffers around to handle correctly.*/
case QR_MODE_BYTE:{
in=(char *)entry->payload.data.buf;
inleft=entry->payload.data.len;
out=sa_text+sa_ntext;
outleft=sa_ctext-sa_ntext;
/*If we have no specified encoding, attempt to auto-detect it.*/
if(eci<0){
int ei;
/*First check for the UTF-8 BOM.*/
if(inleft>=3&&
in[0]==(char)0xEF&&in[1]==(char)0xBB&&in[2]==(char)0xBF){
in+=3;
inleft-=3;
/*Actually try converting (to check validity).*/
err=utf8_cd==(iconv_t)-1||
zbar_iconv((unsigned char* )utf8_cd,&in,&inleft,&out,&outleft)==(size_t)-1;
if(!err){
sa_ntext=out-sa_text;
enc_list_mtf(enc_list,utf8_cd);
continue;
}
in=(char *)entry->payload.data.buf;
inleft=entry->payload.data.len;
out=sa_text+sa_ntext;
outleft=sa_ctext-sa_ntext;
}
/*If the text is 8-bit clean, prefer UTF-8 over SJIS, since SJIS
will corrupt the backslashes used for DoCoMo formats.*/
else if(text_is_ascii((unsigned char *)in,inleft)){
enc_list_mtf(enc_list,utf8_cd);
}
/*Try our list of encodings.*/
for(ei=0;ei<3;ei++)if(enc_list[ei]!=(iconv_t)-1){
/*According to the standard, ISO/IEC 8859-1 (one hyphen) is
supposed to be used, but reality is not always so.
It's got an invalid range that is used often with SJIS
and UTF-8, though, which makes detection easier.
However, iconv() does not properly reject characters in
those ranges, since ISO-8859-1 (two hyphens) defines a
number of seldom-used control code characters there.
So if we see any of those characters, move this
conversion to the end of the list.*/
if(ei<2&&enc_list[ei]==latin1_cd&&
!text_is_latin1((unsigned char *)in,inleft)){
int ej;
for(ej=ei+1;ej<3;ej++)enc_list[ej-1]=enc_list[ej];
enc_list[2]=latin1_cd;
}
err=zbar_iconv((unsigned char* )enc_list[ei],&in,&inleft,&out,&outleft)==(size_t)-1;
if(!err){
sa_ntext=out-sa_text;
enc_list_mtf(enc_list,enc_list[ei]);
break;
}
in=(char *)entry->payload.data.buf;
inleft=entry->payload.data.len;
out=sa_text+sa_ntext;
outleft=sa_ctext-sa_ntext;
}
}
/*We were actually given a character set; use it.*/
else{
err=eci_cd==(iconv_t)-1||
zbar_iconv((unsigned char* )eci_cd,&in,&inleft,&out,&outleft)==(size_t)-1;
if(!err)sa_ntext=out-sa_text;
}
}break;
/*Kanji mode always uses SJIS.*/
case QR_MODE_KANJI:{
in=(char *)entry->payload.data.buf;
inleft=entry->payload.data.len;
out=sa_text+sa_ntext;
outleft=sa_ctext-sa_ntext;
err=sjis_cd==(iconv_t)-1||
zbar_iconv((unsigned char* )sjis_cd,&in,&inleft,&out,&outleft)==(size_t)-1;
if(!err)sa_ntext=out-sa_text;
}break;
/*Check to see if a character set was specified.*/
case QR_MODE_ECI:{
const char *enc;
char buf[16];
unsigned cur_eci;
cur_eci=entry->payload.eci;
if(cur_eci<=QR_ECI_ISO8859_16&&cur_eci!=14){
if(cur_eci!=QR_ECI_GLI0&&cur_eci!=QR_ECI_CP437){
sprintf(buf,"ISO8859-%i",QR_MAXI(cur_eci,3)-2);
enc=buf;
}
/*Note that CP437 requires an iconv compiled with
--enable-extra-encodings, and thus may not be available.*/
else enc="CP437";
}
else if(cur_eci==QR_ECI_SJIS)enc="SJIS";
/*Don't know what this ECI code specifies, but not an encoding that
we recognize.*/
else continue;
eci=cur_eci;
eci_cd=zbar_iconv_open("UTF-8",enc);
}break;
/*Silence stupid compiler warnings.*/
default:break;
}
}
/*If eci should be reset between codes, do so.*/
if(eci<=QR_ECI_GLI1){
eci=-1;
if(eci_cd!=(iconv_t)-1)zbar_iconv_close(eci_cd);
}
}
if(eci_cd!=(iconv_t)-1)zbar_iconv_close(eci_cd);
if(!err){
sa_text[sa_ntext++]='\0';
if(sa_ctext+1>sa_ntext){
sa_text=(char *)realloc(sa_text,sa_ntext*sizeof(*sa_text));
}
zbar_symbol_t *sa_sym;
if(sa_size == 1)
sa_sym = syms;
else {
/* create "virtual" container symbol for composite result */
sa_sym = _zbar_image_scanner_alloc_sym(iscn, ZBAR_QRCODE, 0);
sa_sym->syms = _zbar_symbol_set_create();
sa_sym->syms->head = syms;
/* cheap out w/axis aligned bbox for now */
int xmin = img->width, xmax = -2;
int ymin = img->height, ymax = -2;
/* fixup data references */
for(; syms; syms = syms->next) {
_zbar_symbol_refcnt(syms, 1);
if(syms->type == ZBAR_PARTIAL)
sa_sym->type = ZBAR_PARTIAL;
else
for(j = 0; j < syms->npts; j++) {
int u = syms->pts[j].x;
if(xmin >= u) xmin = u - 1;
if(xmax <= u) xmax = u + 1;
u = syms->pts[j].y;
if(ymin >= u) ymin = u - 1;
if(ymax <= u) ymax = u + 1;
}
syms->data = sa_text + syms->datalen;
int next = (syms->next) ? syms->next->datalen : sa_ntext;
assert(next > syms->datalen);
syms->datalen = next - syms->datalen - 1;
}
if(xmax >= -1) {
sym_add_point(sa_sym, xmin, ymin);
sym_add_point(sa_sym, xmin, ymax);
sym_add_point(sa_sym, xmax, ymax);
sym_add_point(sa_sym, xmax, ymin);
}
}
sa_sym->data = sa_text;
sa_sym->data_alloc = sa_ntext;
sa_sym->datalen = sa_ntext - 1;
_zbar_image_scanner_add_sym(iscn, sa_sym);
}
else {
_zbar_image_scanner_recycle_syms(iscn, syms);
free(sa_text);
}
}
if(utf8_cd!=(iconv_t)-1)zbar_iconv_close(utf8_cd);
if(sjis_cd!=(iconv_t)-1)zbar_iconv_close(sjis_cd);
if(latin1_cd!=(iconv_t)-1)zbar_iconv_close(latin1_cd);
free(mark);
return ntext;
}
+799
View File
@@ -0,0 +1,799 @@
/*Copyright (C) 1991-1995 Henry Minsky (hqm@ua.com, hqm@ai.mit.edu)
Copyright (C) 2008-2009 Timothy B. Terriberry (tterribe@xiph.org)
You can redistribute this library and/or modify it under the terms of the
GNU Lesser General Public License as published by the Free Software
Foundation; either version 2.1 of the License, or (at your option) any later
version.*/
#include <stdlib.h>
#include <string.h>
#include "rs.h"
/*Reed-Solomon encoder and decoder.
Original implementation (C) Henry Minsky (hqm@ua.com, hqm@ai.mit.edu),
Universal Access 1991-1995.
Updates by Timothy B. Terriberry (C) 2008-2009:
- Properly reject codes when error-locator polynomial has repeated roots or
non-trivial irreducible factors.
- Removed the hard-coded parity size and performed general API cleanup.
- Allow multiple representations of GF(2**8), since different standards use
different irreducible polynomials.
- Allow different starting indices for the generator polynomial, since
different standards use different values.
- Greatly reduced the computation by eliminating unnecessary operations.
- Explicitly solve for the roots of low-degree polynomials instead of using
an exhaustive search.
This is another major speed boost when there are few errors.*/
/*Galois Field arithmetic in GF(2**8).*/
void rs_gf256_init(rs_gf256 *_gf,unsigned _ppoly){
unsigned p;
int i;
/*Initialize the table of powers of a primtive root, alpha=0x02.*/
p=1;
for(i=0;i<256;i++){
_gf->exp[i]=_gf->exp[i+255]=p;
p=((p<<1)^(-(p>>7)&_ppoly))&0xFF;
}
/*Invert the table to recover the logs.*/
for(i=0;i<255;i++)_gf->log[_gf->exp[i]]=i;
/*Note that we rely on the fact that _gf->log[0]=0 below.*/
_gf->log[0]=0;
}
/*Multiplication in GF(2**8) using logarithms.*/
static unsigned rs_gmul(const rs_gf256 *_gf,unsigned _a,unsigned _b){
return _a==0||_b==0?0:_gf->exp[_gf->log[_a]+_gf->log[_b]];
}
/*Division in GF(2**8) using logarithms.
The result of division by zero is undefined.*/
static unsigned rs_gdiv(const rs_gf256 *_gf,unsigned _a,unsigned _b){
return _a==0?0:_gf->exp[_gf->log[_a]+255-_gf->log[_b]];
}
/*Multiplication in GF(2**8) when one of the numbers is known to be non-zero
(proven by representing it by its logarithm).*/
static unsigned rs_hgmul(const rs_gf256 *_gf,unsigned _a,unsigned _logb){
return _a==0?0:_gf->exp[_gf->log[_a]+_logb];
}
/*Square root in GF(2**8) using logarithms.*/
static unsigned rs_gsqrt(const rs_gf256 *_gf,unsigned _a){
unsigned loga;
if(!_a)return 0;
loga=_gf->log[_a];
return _gf->exp[(loga+(255&-(loga&1)))>>1];
}
/*Polynomial root finding in GF(2**8).
Each routine returns a list of the distinct roots (i.e., with duplicates
removed).*/
/*Solve a quadratic equation x**2 + _b*x + _c in GF(2**8) using the method
of~\cite{Wal99}.
Returns the number of distinct roots.
ARTICLE{Wal99,
author="C. Wayne Walker",
title="New Formulas for Solving Quadratic Equations over Certain Finite
Fields",
journal="{IEEE} Transactions on Information Theory",
volume=45,
number=1,
pages="283--284",
month=Jan,
year=1999
}*/
static int rs_quadratic_solve(const rs_gf256 *_gf,unsigned _b,unsigned _c,
unsigned char _x[2]){
unsigned b;
unsigned logb;
unsigned logb2;
unsigned logb4;
unsigned logb8;
unsigned logb12;
unsigned logb14;
unsigned logc;
unsigned logc2;
unsigned logc4;
unsigned c8;
unsigned g3;
unsigned z3;
unsigned l3;
unsigned c0;
unsigned g2;
unsigned l2;
unsigned z2;
int inc;
/*If _b is zero, all we need is a square root.*/
if(!_b){
_x[0]=rs_gsqrt(_gf,_c);
return 1;
}
/*If _c is zero, 0 and _b are the roots.*/
if(!_c){
_x[0]=0;
_x[1]=_b;
return 2;
}
logb=_gf->log[_b];
logc=_gf->log[_c];
/*If _b lies in GF(2**4), scale x to move it out.*/
inc=logb%(255/15)==0;
if(inc){
b=_gf->exp[logb+254];
logb=_gf->log[b];
_c=_gf->exp[logc+253];
logc=_gf->log[_c];
}
else b=_b;
logb2=_gf->log[_gf->exp[logb<<1]];
logb4=_gf->log[_gf->exp[logb2<<1]];
logb8=_gf->log[_gf->exp[logb4<<1]];
logb12=_gf->log[_gf->exp[logb4+logb8]];
logb14=_gf->log[_gf->exp[logb2+logb12]];
logc2=_gf->log[_gf->exp[logc<<1]];
logc4=_gf->log[_gf->exp[logc2<<1]];
c8=_gf->exp[logc4<<1];
g3=rs_hgmul(_gf,
_gf->exp[logb14+logc]^_gf->exp[logb12+logc2]^_gf->exp[logb8+logc4]^c8,logb);
/*If g3 doesn't lie in GF(2**4), then our roots lie in an extension field.
Note that we rely on the fact that _gf->log[0]==0 here.*/
if(_gf->log[g3]%(255/15)!=0)return 0;
/*Construct the corresponding quadratic in GF(2**4):
x**2 + x/alpha**(255/15) + l3/alpha**(2*(255/15))*/
z3=rs_gdiv(_gf,g3,_gf->exp[logb8<<1]^b);
l3=rs_hgmul(_gf,rs_gmul(_gf,z3,z3)^rs_hgmul(_gf,z3,logb)^_c,255-logb2);
c0=rs_hgmul(_gf,l3,255-2*(255/15));
/*Construct the corresponding quadratic in GF(2**2):
x**2 + x/alpha**(255/3) + l2/alpha**(2*(255/3))*/
g2=rs_hgmul(_gf,
rs_hgmul(_gf,c0,255-2*(255/15))^rs_gmul(_gf,c0,c0),255-255/15);
z2=rs_gdiv(_gf,g2,_gf->exp[255-(255/15)*4]^_gf->exp[255-(255/15)]);
l2=rs_hgmul(_gf,
rs_gmul(_gf,z2,z2)^rs_hgmul(_gf,z2,255-(255/15))^c0,2*(255/15));
/*Back substitute to the solution in the original field.*/
_x[0]=_gf->exp[_gf->log[z3^rs_hgmul(_gf,
rs_hgmul(_gf,l2,255/3)^rs_hgmul(_gf,z2,255/15),logb)]+inc];
_x[1]=_x[0]^_b;
return 2;
}
/*Solve a cubic equation x**3 + _a*x**2 + _b*x + _c in GF(2**8).
Returns the number of distinct roots.*/
static int rs_cubic_solve(const rs_gf256 *_gf,
unsigned _a,unsigned _b,unsigned _c,unsigned char _x[3]){
unsigned k;
unsigned logd;
unsigned d2;
unsigned logd2;
unsigned logw;
int nroots;
/*If _c is zero, factor out the 0 root.*/
if(!_c){
nroots=rs_quadratic_solve(_gf,_a,_b,_x);
if(_b)_x[nroots++]=0;
return nroots;
}
/*Substitute x=_a+y*sqrt(_a**2+_b) to get y**3 + y + k == 0,
k = (_a*_b+c)/(_a**2+b)**(3/2).*/
k=rs_gmul(_gf,_a,_b)^_c;
d2=rs_gmul(_gf,_a,_a)^_b;
if(!d2){
int logx;
if(!k){
/*We have a triple root.*/
_x[0]=_a;
return 1;
}
logx=_gf->log[k];
if(logx%3!=0)return 0;
logx/=3;
_x[0]=_a^_gf->exp[logx];
_x[1]=_a^_gf->exp[logx+255/3];
_x[2]=_a^_x[0]^_x[1];
return 3;
}
logd2=_gf->log[d2];
logd=(logd2+(255&-(logd2&1)))>>1;
k=rs_gdiv(_gf,k,_gf->exp[logd+logd2]);
/*Substitute y=w+1/w and z=w**3 to get z**2 + k*z + 1 == 0.*/
nroots=rs_quadratic_solve(_gf,k,1,_x);
if(nroots<1){
/*The Reed-Solomon code is only valid if we can find 3 distinct roots in
GF(2**8), so if we know there's only one, we don't actually need to find
it.
Note that we're also called by the quartic solver, but if we contain a
non-trivial irreducible factor, than so does the original
quartic~\cite{LW72}, and failing to return a root here actually saves us
some work there, also.*/
return 0;
}
/*Recover w from z.*/
logw=_gf->log[_x[0]];
if(logw){
if(logw%3!=0)return 0;
logw/=3;
/*Recover x from w.*/
_x[0]=_gf->exp[_gf->log[_gf->exp[logw]^_gf->exp[255-logw]]+logd]^_a;
logw+=255/3;
_x[1]=_gf->exp[_gf->log[_gf->exp[logw]^_gf->exp[255-logw]]+logd]^_a;
_x[2]=_x[0]^_x[1]^_a;
return 3;
}
else{
_x[0]=_a;
/*In this case _x[1] is a double root, so we know the Reed-Solomon code is
invalid.
Note that we still have to return at least one root, because if we're
being called by the quartic solver, the quartic might still have 4
distinct roots.
But we don't need more than one root, so we can avoid computing the
expensive one.*/
/*_x[1]=_gf->exp[_gf->log[_gf->exp[255/3]^_gf->exp[2*(255/3)]]+logd]^_a;*/
return 1;
}
}
/*Solve a quartic equation x**4 + _a*x**3 + _b*x**2 + _c*x + _d in GF(2**8) by
decomposing it into the cases given by~\cite{LW72}.
Returns the number of distinct roots.
@ARTICLE{LW72,
author="Philip A. Leonard and Kenneth S. Williams",
title="Quartics over $GF(2^n)$",
journal="Proceedings of the American Mathematical Society",
volume=36,
number=2,
pages="347--450",
month=Dec,
year=1972
}*/
static int rs_quartic_solve(const rs_gf256 *_gf,
unsigned _a,unsigned _b,unsigned _c,unsigned _d,unsigned char _x[3]){
unsigned r;
unsigned s;
unsigned t;
unsigned b;
int nroots;
int i;
/*If _d is zero, factor out the 0 root.*/
if(!_d){
nroots=rs_cubic_solve(_gf,_a,_b,_c,_x);
if(_c)_x[nroots++]=0;
return nroots;
}
if(_a){
unsigned loga;
/*Substitute x=(1/y) + sqrt(_c/_a) to eliminate the cubic term.*/
loga=_gf->log[_a];
r=rs_hgmul(_gf,_c,255-loga);
s=rs_gsqrt(_gf,r);
t=_d^rs_gmul(_gf,_b,r)^rs_gmul(_gf,r,r);
if(t){
unsigned logti;
logti=255-_gf->log[t];
/*The result is still quartic, but with no cubic term.*/
nroots=rs_quartic_solve(_gf,0,rs_hgmul(_gf,_b^rs_hgmul(_gf,s,loga),logti),
_gf->exp[loga+logti],_gf->exp[logti],_x);
for(i=0;i<nroots;i++)_x[i]=_gf->exp[255-_gf->log[_x[i]]]^s;
}
else{
/*s must be a root~\cite{LW72}, and is in fact a double-root~\cite{CCO69}.
Thus we're left with only a quadratic to solve.
@ARTICLE{CCO69,
author="Robert T. Chien and B. D. Cunningham and I. B. Oldham",
title="Hybrid Methods for Finding Roots of a Polynomial---With
Applications to {BCH} Decoding",
journal="{IEEE} Transactions on Information Theory",
volume=15,
number=2,
pages="329--335",
month=Mar,
year=1969
}*/
nroots=rs_quadratic_solve(_gf,_a,_b^r,_x);
/*s may be a triple root if s=_b/_a, but not quadruple, since _a!=0.*/
if(nroots!=2||(_x[0]!=s&&_x[1]!=s))_x[nroots++]=s;
}
return nroots;
}
/*If there are no odd powers, it's really just a quadratic in disguise.*/
if(!_c)return rs_quadratic_solve(_gf,rs_gsqrt(_gf,_b),rs_gsqrt(_gf,_d),_x);
/*Factor into (x**2 + r*x + s)*(x**2 + r*x + t) by solving for r, which can
be shown to satisfy r**3 + _b*r + _c == 0.*/
nroots=rs_cubic_solve(_gf,0,_b,_c,_x);
if(nroots<1){
/*The Reed-Solomon code is only valid if we can find 4 distinct roots in
GF(2**8).
If the cubic does not factor into 3 (possibly duplicate) roots, then we
know that the quartic must have a non-trivial irreducible factor.*/
return 0;
}
r=_x[0];
/*Now solve for s and t.*/
b=rs_gdiv(_gf,_c,r);
nroots=rs_quadratic_solve(_gf,b,_d,_x);
if(nroots<2)return 0;
s=_x[0];
t=_x[1];
/*_c=r*(s^t) was non-zero, so s and t must be distinct.
But if z is a root of z**2 ^ r*z ^ s, then so is (z^r), and s = z*(z^r).
Hence if z is also a root of z**2 ^ r*z ^ t, then t = s, a contradiction.
Thus all four roots are distinct, if they exist.*/
nroots=rs_quadratic_solve(_gf,r,s,_x);
return nroots+rs_quadratic_solve(_gf,r,t,_x+nroots);
}
/*Polynomial arithmetic with coefficients in GF(2**8).*/
static void rs_poly_zero(unsigned char *_p,int _dp1){
memset(_p,0,_dp1*sizeof(*_p));
}
static void rs_poly_copy(unsigned char *_p,const unsigned char *_q,int _dp1){
memcpy(_p,_q,_dp1*sizeof(*_p));
}
/*Multiply the polynomial by the free variable, x (shift the coefficients).
The number of coefficients, _dp1, must be non-zero.*/
static void rs_poly_mul_x(unsigned char *_p,const unsigned char *_q,int _dp1){
memmove(_p+1,_q,(_dp1-1)*sizeof(*_p));
_p[0]=0;
}
/*Divide the polynomial by the free variable, x (shift the coefficients).
The number of coefficients, _dp1, must be non-zero.*/
static void rs_poly_div_x(unsigned char *_p,const unsigned char *_q,int _dp1){
memmove(_p,_q+1,(_dp1-1)*sizeof(*_p));
_p[_dp1-1]=0;
}
/*Compute the first (d+1) coefficients of the product of a degree e and a
degree f polynomial.*/
static void rs_poly_mult(const rs_gf256 *_gf,unsigned char *_p,int _dp1,
const unsigned char *_q,int _ep1,const unsigned char *_r,int _fp1){
int m;
int i;
rs_poly_zero(_p,_dp1);
m=_ep1<_dp1?_ep1:_dp1;
for(i=0;i<m;i++)if(_q[i]!=0){
unsigned logqi;
int n;
int j;
n=_dp1-i<_fp1?_dp1-i:_fp1;
logqi=_gf->log[_q[i]];
for(j=0;j<n;j++)_p[i+j]^=rs_hgmul(_gf,_r[j],logqi);
}
}
/*Decoding.*/
/*Computes the syndrome of a codeword.*/
static void rs_calc_syndrome(const rs_gf256 *_gf,int _m0,
unsigned char *_s,int _npar,const unsigned char *_data,int _ndata){
int i;
int j;
for(j=0;j<_npar;j++){
unsigned alphaj;
unsigned sj;
sj=0;
alphaj=_gf->log[_gf->exp[j+_m0]];
for(i=0;i<_ndata;i++)sj=_data[i]^rs_hgmul(_gf,sj,alphaj);
_s[j]=sj;
}
}
/*Berlekamp-Peterson and Berlekamp-Massey Algorithms for error-location,
modified to handle known erasures, from \cite{CC81}, p. 205.
This finds the coefficients of the error locator polynomial.
The roots are then found by looking for the values of alpha**n where
evaluating the polynomial yields zero.
Error correction is done using the error-evaluator equation on p. 207.
@BOOK{CC81,
author="George C. Clark, Jr and J. Bibb Cain",
title="Error-Correction Coding for Digitial Communications",
series="Applications of Communications Theory",
publisher="Springer",
address="New York, NY",
month=Jun,
year=1981
}*/
/*Initialize lambda to the product of (1-x*alpha**e[i]) for erasure locations
e[i].
Note that the user passes in array indices counting from the beginning of the
data, while our polynomial indexes starting from the end, so
e[i]=(_ndata-1)-_erasures[i].*/
static void rs_init_lambda(const rs_gf256 *_gf,unsigned char *_lambda,int _npar,
const unsigned char *_erasures,int _nerasures,int _ndata){
int i;
int j;
rs_poly_zero(_lambda,(_npar<4?4:_npar)+1);
_lambda[0]=1;
for(i=0;i<_nerasures;i++)for(j=i+1;j>0;j--){
_lambda[j]^=rs_hgmul(_gf,_lambda[j-1],_ndata-1-_erasures[i]);
}
}
/*From \cite{CC81}, p. 216.
Returns the number of errors detected (degree of _lambda).*/
static int rs_modified_berlekamp_massey(const rs_gf256 *_gf,
unsigned char *_lambda,const unsigned char *_s,unsigned char *_omega,int _npar,
const unsigned char *_erasures,int _nerasures,int _ndata){
unsigned char tt[256];
int n;
int l;
int k;
int i;
/*Initialize _lambda, the error locator-polynomial, with the location of
known erasures.*/
rs_init_lambda(_gf,_lambda,_npar,_erasures,_nerasures,_ndata);
rs_poly_copy(tt,_lambda,_npar+1);
l=_nerasures;
k=0;
for(n=_nerasures+1;n<=_npar;n++){
unsigned d;
rs_poly_mul_x(tt,tt,n-k+1);
d=0;
for(i=0;i<=l;i++)d^=rs_gmul(_gf,_lambda[i],_s[n-1-i]);
if(d!=0){
unsigned logd;
logd=_gf->log[d];
if(l<n-k){
int t;
for(i=0;i<=n-k;i++){
unsigned tti;
tti=tt[i];
tt[i]=rs_hgmul(_gf,_lambda[i],255-logd);
_lambda[i]=_lambda[i]^rs_hgmul(_gf,tti,logd);
}
t=n-k;
k=n-l;
l=t;
}
else for(i=0;i<=l;i++)_lambda[i]=_lambda[i]^rs_hgmul(_gf,tt[i],logd);
}
}
rs_poly_mult(_gf,_omega,_npar,_lambda,l+1,_s,_npar);
return l;
}
/*Finds all the roots of an error-locator polynomial _lambda by evaluating it
at successive values of alpha, and returns the positions of the associated
errors in _epos.
Returns the number of valid roots identified.*/
static int rs_find_roots(const rs_gf256 *_gf,unsigned char *_epos,
const unsigned char *_lambda,int _nerrors,int _ndata){
unsigned alpha;
int nroots;
int i;
nroots=0;
if(_nerrors<=4){
/*Explicit solutions for higher degrees are possible.
Chien uses large lookup tables to solve quintics, and Truong et al. give
special algorithms for degree up through 11, though they use exhaustive
search (with reduced complexity) for some portions.
Quartics are good enough for reading CDs, and represent a reasonable code
complexity trade-off without requiring any extra tables.
Note that _lambda[0] is always 1.*/
_nerrors=rs_quartic_solve(_gf,_lambda[1],_lambda[2],_lambda[3],_lambda[4],
_epos);
for(i=0;i<_nerrors;i++)if(_epos[i]){
alpha=_gf->log[_epos[i]];
if((int)alpha<_ndata)_epos[nroots++]=alpha;
}
return nroots;
}
else for(alpha=0;(int)alpha<_ndata;alpha++){
unsigned alphai;
unsigned sum;
sum=0;
alphai=0;
for(i=0;i<=_nerrors;i++){
sum^=rs_hgmul(_gf,_lambda[_nerrors-i],alphai);
alphai=_gf->log[_gf->exp[alphai+alpha]];
}
if(!sum)_epos[nroots++]=alpha;
}
return nroots;
}
/*Corrects a codeword with _ndata<256 bytes, of which the last _npar are parity
bytes.
Known locations of errors can be passed in the _erasures array.
Twice as many (up to _npar) errors with a known location can be corrected
compared to errors with an unknown location.
Returns the number of errors corrected if successful, or a negative number if
the message could not be corrected because too many errors were detected.*/
int rs_correct(const rs_gf256 *_gf,int _m0,unsigned char *_data,int _ndata,
int _npar,const unsigned char *_erasures,int _nerasures){
/*lambda must have storage for at least five entries to avoid special cases
in the low-degree polynomial solver.*/
unsigned char lambda[256];
unsigned char omega[256];
unsigned char epos[256];
unsigned char s[256];
int i;
/*If we already have too many erasures, we can't possibly succeed.*/
if(_nerasures>_npar)return -1;
/*Compute the syndrome values.*/
rs_calc_syndrome(_gf,_m0,s,_npar,_data,_ndata);
/*Check for a non-zero value.*/
for(i=0;i<_npar;i++)if(s[i]){
int nerrors;
int j;
/*Construct the error locator polynomial.*/
nerrors=rs_modified_berlekamp_massey(_gf,lambda,s,omega,_npar,
_erasures,_nerasures,_ndata);
/*If we can't locate any errors, we can't force the syndrome values to
zero, and must have a decoding error.
Conversely, if we have too many errors, there's no reason to even attempt
the root search.*/
if(nerrors<=0||nerrors-_nerasures>(_npar-_nerasures)>>1)return -1;
/*Compute the locations of the errors.
If they are not all distinct, or some of them were outside the valid
range for our block size, we have a decoding error.*/
if(rs_find_roots(_gf,epos,lambda,nerrors,_ndata)<nerrors)return -1;
/*Now compute the error magnitudes.*/
for(i=0;i<nerrors;i++){
unsigned a;
unsigned b;
unsigned alpha;
unsigned alphan1;
unsigned alphan2;
unsigned alphanj;
alpha=epos[i];
/*Evaluate omega at alpha**-1.*/
a=0;
alphan1=255-alpha;
alphanj=0;
for(j=0;j<_npar;j++){
a^=rs_hgmul(_gf,omega[j],alphanj);
alphanj=_gf->log[_gf->exp[alphanj+alphan1]];
}
/*Evaluate the derivative of lambda at alpha**-1
All the odd powers vanish.*/
b=0;
alphan2=_gf->log[_gf->exp[alphan1<<1]];
alphanj=alphan1+_m0*alpha%255;
for(j=1;j<=_npar;j+=2){
b^=rs_hgmul(_gf,lambda[j],alphanj);
alphanj=_gf->log[_gf->exp[alphanj+alphan2]];
}
/*Apply the correction.*/
_data[_ndata-1-alpha]^=rs_gdiv(_gf,a,b);
}
return nerrors;
}
return 0;
}
/*Encoding.*/
/*Create an _npar-coefficient generator polynomial for a Reed-Solomon code
with _npar<256 parity bytes.*/
void rs_compute_genpoly(const rs_gf256 *_gf,int _m0,
unsigned char *_genpoly,int _npar){
int i;
if(_npar<=0)return;
rs_poly_zero(_genpoly,_npar);
_genpoly[0]=1;
/*Multiply by (x+alpha^i) for i = 1 ... _ndata.*/
for(i=0;i<_npar;i++){
unsigned alphai;
int n;
int j;
n=i+1<_npar-1?i+1:_npar-1;
alphai=_gf->log[_gf->exp[_m0+i]];
for(j=n;j>0;j--)_genpoly[j]=_genpoly[j-1]^rs_hgmul(_gf,_genpoly[j],alphai);
_genpoly[0]=rs_hgmul(_gf,_genpoly[0],alphai);
}
}
/*Adds _npar<=_ndata parity bytes to an _ndata-_npar byte message.
_data must contain room for _ndata<256 bytes.*/
void rs_encode(const rs_gf256 *_gf,unsigned char *_data,int _ndata,
const unsigned char *_genpoly,int _npar){
unsigned char *lfsr;
unsigned d;
int i;
int j;
if(_npar<=0)return;
lfsr=_data+_ndata-_npar;
rs_poly_zero(lfsr,_npar);
for(i=0;i<_ndata-_npar;i++){
d=_data[i]^lfsr[0];
if(d){
unsigned logd;
logd=_gf->log[d];
for(j=0;j<_npar-1;j++){
lfsr[j]=lfsr[j+1]^rs_hgmul(_gf,_genpoly[_npar-1-j],logd);
}
lfsr[_npar-1]=rs_hgmul(_gf,_genpoly[0],logd);
}
else rs_poly_div_x(lfsr,lfsr,_npar);
}
}
#if defined(RS_TEST_ENC)
#include <stdio.h>
#include <stdlib.h>
int main(void){
rs_gf256 gf;
int k;
rs_gf256_init(&gf,QR_PPOLY);
srand(0);
for(k=0;k<64*1024;k++){
unsigned char genpoly[256];
unsigned char data[256];
unsigned char epos[256];
int ndata;
int npar;
int nerrors;
int i;
ndata=rand()&0xFF;
npar=ndata>0?rand()%ndata:0;
for(i=0;i<ndata-npar;i++)data[i]=rand()&0xFF;
rs_compute_genpoly(&gf,QR_M0,genpoly,npar);
rs_encode(&gf,QR_M0,data,ndata,genpoly,npar);
/*Write a clean version of the codeword.*/
printf("%i %i",ndata,npar);
for(i=0;i<ndata;i++)printf(" %i",data[i]);
printf(" 0\n");
/*Write the correct output to compare the decoder against.*/
fprintf(stderr,"Success!\n",nerrors);
for(i=0;i<ndata;i++)fprintf(stderr,"%i%s",data[i],i+1<ndata?" ":"\n");
if(npar>0){
/*Corrupt it.*/
nerrors=rand()%(npar+1);
if(nerrors>0){
/*This test is not quite correct: there could be so many errors it
comes within (npar>>1) errors of another valid codeword.
I don't know a simple way to test for that without trying to decode
the corrupt codeword, though, which is the very code we're trying to
test.*/
if(nerrors<=npar>>1){
fprintf(stderr,"Success!\n",nerrors);
for(i=0;i<ndata;i++){
fprintf(stderr,"%i%s",data[i],i+1<ndata?" ":"\n");
}
}
else fprintf(stderr,"Failure.\n");
fprintf(stderr,"Success!\n",nerrors);
for(i=0;i<ndata;i++)fprintf(stderr,"%i%s",data[i],i+1<ndata?" ":"\n");
for(i=0;i<ndata;i++)epos[i]=i;
for(i=0;i<nerrors;i++){
unsigned char e;
int ei;
ei=rand()%(ndata-i)+i;
e=epos[ei];
epos[ei]=epos[i];
epos[i]=e;
data[e]^=rand()%255+1;
}
/*First with no erasure locations.*/
printf("%i %i",ndata,npar);
for(i=0;i<ndata;i++)printf(" %i",data[i]);
printf(" 0\n");
/*Now with erasure locations.*/
printf("%i %i",ndata,npar);
for(i=0;i<ndata;i++)printf(" %i",data[i]);
printf(" %i",nerrors);
for(i=0;i<nerrors;i++)printf(" %i",epos[i]);
printf("\n");
}
}
}
return 0;
}
#endif
#if defined(RS_TEST_DEC)
#include <stdio.h>
int main(void){
rs_gf256 gf;
rs_gf256_init(&gf,QR_PPOLY);
for(;;){
unsigned char data[255];
unsigned char erasures[255];
int idata[255];
int ierasures[255];
int ndata;
int npar;
int nerasures;
int nerrors;
int i;
if(scanf("%i",&ndata)<1||ndata<0||ndata>255||
scanf("%i",&npar)<1||npar<0||npar>ndata)break;
for(i=0;i<ndata;i++){
if(scanf("%i",idata+i)<1||idata[i]<0||idata[i]>255)break;
data[i]=idata[i];
}
if(i<ndata)break;
if(scanf("%i",&nerasures)<1||nerasures<0||nerasures>ndata)break;
for(i=0;i<nerasures;i++){
if(scanf("%i",ierasures+i)<1||ierasures[i]<0||ierasures[i]>=ndata)break;
erasures[i]=ierasures[i];
}
nerrors=rs_correct(&gf,QR_M0,data,ndata,npar,erasures,nerasures);
if(nerrors>=0){
unsigned char data2[255];
unsigned char genpoly[255];
for(i=0;i<ndata-npar;i++)data2[i]=data[i];
rs_compute_genpoly(&gf,QR_M0,genpoly,npar);
rs_encode(&gf,QR_M0,data2,ndata,genpoly,npar);
for(i=ndata-npar;i<ndata;i++)if(data[i]!=data2[i]){
printf("Abject failure! %i!=%i\n",data[i],data2[i]);
}
printf("Success!\n",nerrors);
for(i=0;i<ndata;i++)printf("%i%s",data[i],i+1<ndata?" ":"\n");
}
else printf("Failure.\n");
}
return 0;
}
#endif
#if defined(RS_TEST_ROOTS)
#include <stdio.h>
/*Exhaustively test the root finder.*/
int main(void){
rs_gf256 gf;
int a;
int b;
int c;
int d;
rs_gf256_init(&gf,QR_PPOLY);
for(a=0;a<256;a++)for(b=0;b<256;b++)for(c=0;c<256;c++)for(d=0;d<256;d++){
unsigned char x[4];
unsigned char r[4];
unsigned x2;
unsigned e[5];
int nroots;
int mroots;
int i;
int j;
nroots=rs_quartic_solve(&gf,a,b,c,d,x);
for(i=0;i<nroots;i++){
x2=rs_gmul(&gf,x[i],x[i]);
e[0]=rs_gmul(&gf,x2,x2)^rs_gmul(&gf,a,rs_gmul(&gf,x[i],x2))^
rs_gmul(&gf,b,x2)^rs_gmul(&gf,c,x[i])^d;
if(e[0]){
printf("Invalid root: (0x%02X)**4 ^ 0x%02X*(0x%02X)**3 ^ "
"0x%02X*(0x%02X)**2 ^ 0x%02X(0x%02X) ^ 0x%02X = 0x%02X\n",
x[i],a,x[i],b,x[i],c,x[i],d,e[0]);
}
for(j=0;j<i;j++)if(x[i]==x[j]){
printf("Repeated root %i=%i: (0x%02X)**4 ^ 0x%02X*(0x%02X)**3 ^ "
"0x%02X*(0x%02X)**2 ^ 0x%02X(0x%02X) ^ 0x%02X = 0x%02X\n",
i,j,x[i],a,x[i],b,x[i],c,x[i],d,e[0]);
}
}
mroots=0;
for(j=1;j<256;j++){
int logx;
int logx2;
logx=gf.log[j];
logx2=gf.log[gf.exp[logx<<1]];
e[mroots]=d^rs_hgmul(&gf,c,logx)^rs_hgmul(&gf,b,logx2)^
rs_hgmul(&gf,a,gf.log[gf.exp[logx+logx2]])^gf.exp[logx2<<1];
if(!e[mroots])r[mroots++]=j;
}
/*We only care about missing roots if the quartic had 4 distinct, non-zero
roots.*/
if(mroots==4)for(j=0;j<mroots;j++){
for(i=0;i<nroots;i++)if(x[i]==r[j])break;
if(i>=nroots){
printf("Missing root: (0x%02X)**4 ^ 0x%02X*(0x%02X)**3 ^ "
"0x%02X*(0x%02X)**2 ^ 0x%02X(0x%02X) ^ 0x%02X = 0x%02X\n",
r[j],a,r[j],b,r[j],c,r[j],d,e[j]);
}
}
}
return 0;
}
#endif
+66
View File
@@ -0,0 +1,66 @@
/*Copyright (C) 1991-1995 Henry Minsky (hqm@ua.com, hqm@ai.mit.edu)
Copyright (C) 2008-2009 Timothy B. Terriberry (tterribe@xiph.org)
You can redistribute this library and/or modify it under the terms of the
GNU Lesser General Public License as published by the Free Software
Foundation; either version 2.1 of the License, or (at your option) any later
version.*/
#if !defined(_qrcode_rs_H)
# define _qrcode_rs_H (1)
/*This is one of 16 irreducible primitive polynomials of degree 8:
x**8+x**4+x**3+x**2+1.
Under such a polynomial, x (i.e., 0x02) is a generator of GF(2**8).
The high order 1 bit is implicit.
From~\cite{MD88}, Ch. 5, p. 275 by Patel.
@BOOK{MD88,
author="C. Dennis Mee and Eric D. Daniel",
title="Video, Audio, and Instrumentation Recording",
series="Magnetic Recording",
volume=3,
publisher="McGraw-Hill Education",
address="Columbus, OH",
month=Jun,
year=1988
}*/
#define QR_PPOLY (0x1D)
/*The index to start the generator polynomial from (0...254).*/
#define QR_M0 (0)
typedef struct rs_gf256 rs_gf256;
struct rs_gf256{
/*A logarithm table in GF(2**8).*/
unsigned char log[256];
/*An exponential table in GF(2**8): exp[i] contains x^i reduced modulo the
irreducible primitive polynomial used to define the field.
The extra 256 entries are used to do arithmetic mod 255, since some extra
table lookups are generally faster than doing the modulus.*/
unsigned char exp[511];
};
/*Initialize discrete logarithm tables for GF(2**8) using a given primitive
irreducible polynomial.*/
void rs_gf256_init(rs_gf256 *_gf,unsigned _ppoly);
/*Corrects a codeword with _ndata<256 bytes, of which the last _npar are parity
bytes.
Known locations of errors can be passed in the _erasures array.
Twice as many (up to _npar) errors with a known location can be corrected
compared to errors with an unknown location.
Returns the number of errors corrected if successful, or a negative number if
the message could not be corrected because too many errors were detected.*/
int rs_correct(const rs_gf256 *_gf,int _m0,unsigned char *_data,int _ndata,
int _npar,const unsigned char *_erasures,int _nerasures);
/*Create an _npar-coefficient generator polynomial for a Reed-Solomon code with
_npar<256 parity bytes.*/
void rs_compute_genpoly(const rs_gf256 *_gf,int _m0,
unsigned char *_genpoly,int _npar);
/*Adds _npar<=_ndata parity bytes to an _ndata-_npar byte message.
_data must contain room for _ndata<256 bytes.*/
void rs_encode(const rs_gf256 *_gf,unsigned char *_data,int _ndata,
const unsigned char *_genpoly,int _npar);
#endif
+140
View File
@@ -0,0 +1,140 @@
/*Copyright (C) 2008-2009 Timothy B. Terriberry (tterribe@xiph.org)
You can redistribute this library and/or modify it under the terms of the
GNU Lesser General Public License as published by the Free Software
Foundation; either version 2.1 of the License, or (at your option) any later
version.*/
#include <stdlib.h>
#include "util.h"
/*Computes floor(sqrt(_val)) exactly.*/
unsigned qr_isqrt(unsigned _val){
unsigned g;
unsigned b;
int bshift;
/*Uses the second method from
http://www.azillionmonkeys.com/qed/sqroot.html
The main idea is to search for the largest binary digit b such that
(g+b)*(g+b) <= _val, and add it to the solution g.*/
g=0;
b=0x8000;
for(bshift=16;bshift-->0;){
unsigned t;
t=((g<<1)+b)<<bshift;
if(t<=_val){
g+=b;
_val-=t;
}
b>>=1;
}
return g;
}
/*Computes sqrt(_x*_x+_y*_y) using CORDIC.
This implementation is valid for all 32-bit inputs and returns a result
accurate to about 27 bits of precision.
It has been tested for all postiive 16-bit inputs, where it returns correctly
rounded results in 99.998% of cases and the maximum error is
0.500137134862598032 (for _x=48140, _y=63018).
Very nearly all results less than (1<<16) are correctly rounded.
All Pythagorean triples with a hypotenuse of less than ((1<<27)-1) evaluate
correctly, and the total bias over all Pythagorean triples is -0.04579, with
a relative RMS error of 7.2864E-10 and a relative peak error of 7.4387E-9.*/
unsigned qr_ihypot(int _x,int _y){
unsigned x;
unsigned y;
int mask;
int shift;
int u;
int v;
int i;
x=_x=abs(_x);
y=_y=abs(_y);
mask=-(x>y)&(_x^_y);
x^=mask;
y^=mask;
_y^=mask;
shift=31-qr_ilog(y);
shift=QR_MAXI(shift,0);
x=(unsigned)((x<<shift)*0x9B74EDAAULL>>32);
_y=(int)((_y<<shift)*0x9B74EDA9LL>>32);
u=x;
mask=-(_y<0);
x+=(_y+mask)^mask;
_y-=(u+mask)^mask;
u=(x+1)>>1;
v=(_y+1)>>1;
mask=-(_y<0);
x+=(v+mask)^mask;
_y-=(u+mask)^mask;
for(i=1;i<16;i++){
int r;
u=(x+1)>>2;
r=(1<<2*i)>>1;
v=(_y+r)>>2*i;
mask=-(_y<0);
x+=(v+mask)^mask;
_y=(_y-((u+mask)^mask))<<1;
}
return (x+((1U<<shift)>>1))>>shift;
}
#if defined(__GNUC__) && defined(HAVE_FEATURES_H)
# include <features.h>
# if __GNUC_PREREQ(3,4)
# include <limits.h>
# if INT_MAX>=2147483647
# define QR_CLZ0 sizeof(unsigned)*CHAR_BIT
# define QR_CLZ(_x) (__builtin_clz(_x))
# elif LONG_MAX>=2147483647L
# define QR_CLZ0 sizeof(unsigned long)*CHAR_BIT
# define QR_CLZ(_x) (__builtin_clzl(_x))
# endif
# endif
#endif
int qr_ilog(unsigned _v){
#if defined(QR_CLZ)
/*Note that __builtin_clz is not defined when _x==0, according to the gcc
documentation (and that of the x86 BSR instruction that implements it), so
we have to special-case it.*/
return QR_CLZ0-QR_CLZ(_v)&-!!_v;
#else
int ret;
int m;
m=!!(_v&0xFFFF0000)<<4;
_v>>=m;
ret=m;
m=!!(_v&0xFF00)<<3;
_v>>=m;
ret|=m;
m=!!(_v&0xF0)<<2;
_v>>=m;
ret|=m;
m=!!(_v&0xC)<<1;
_v>>=m;
ret|=m;
ret|=!!(_v&0x2);
return ret + !!_v;
#endif
}
#if defined(QR_TEST_SQRT)
#include <math.h>
#include <stdio.h>
/*Exhaustively test the integer square root function.*/
int main(void){
unsigned u;
u=0;
do{
unsigned r;
unsigned s;
r=qr_isqrt(u);
s=(int)sqrt(u);
if(r!=s)printf("%u: %u!=%u\n",u,r,s);
u++;
}
while(u);
return 0;
}
#endif
+48
View File
@@ -0,0 +1,48 @@
/*Copyright (C) 2008-2009 Timothy B. Terriberry (tterribe@xiph.org)
You can redistribute this library and/or modify it under the terms of the
GNU Lesser General Public License as published by the Free Software
Foundation; either version 2.1 of the License, or (at your option) any later
version.*/
#if !defined(_qrcode_util_H)
# define _qrcode_util_H (1)
#define QR_MAXI(_a,_b) ((_a)-(((_a)-(_b))&-((_b)>(_a))))
#define QR_MINI(_a,_b) ((_a)+(((_b)-(_a))&-((_b)<(_a))))
#define QR_SIGNI(_x) (((_x)>0)-((_x)<0))
#define QR_SIGNMASK(_x) (-((_x)<0))
/*Unlike copysign(), simply inverts the sign of _a if _b is negative.*/
#define QR_FLIPSIGNI(_a,_b) (((_a)+QR_SIGNMASK(_b))^QR_SIGNMASK(_b))
#define QR_COPYSIGNI(_a,_b) QR_FLIPSIGNI(abs(_a),_b)
/*Divides a signed integer by a positive value with exact rounding.*/
#define QR_DIVROUND(_x,_y) (((_x)+QR_FLIPSIGNI(_y>>1,_x))/(_y))
#define QR_CLAMPI(_a,_b,_c) (QR_MAXI(_a,QR_MINI(_b,_c)))
#define QR_CLAMP255(_x) ((unsigned char)((((_x)<0)-1)&((_x)|-((_x)>255))))
/*Swaps two integers _a and _b if _a>_b.*/
#define QR_SORT2I(_a,_b) \
do{ \
int t__; \
t__=QR_MINI(_a,_b)^(_a); \
(_a)^=t__; \
(_b)^=t__; \
} \
while(0)
#define QR_ILOG0(_v) (!!((_v)&0x2))
#define QR_ILOG1(_v) (((_v)&0xC)?2+QR_ILOG0((_v)>>2):QR_ILOG0(_v))
#define QR_ILOG2(_v) (((_v)&0xF0)?4+QR_ILOG1((_v)>>4):QR_ILOG1(_v))
#define QR_ILOG3(_v) (((_v)&0xFF00)?8+QR_ILOG2((_v)>>8):QR_ILOG2(_v))
#define QR_ILOG4(_v) (((_v)&0xFFFF0000)?16+QR_ILOG3((_v)>>16):QR_ILOG3(_v))
/*Computes the integer logarithm of a (positive, 32-bit) constant.*/
#define QR_ILOG(_v) ((int)QR_ILOG4((unsigned)(_v)))
/*Multiplies 32-bit numbers _a and _b, adds (possibly 64-bit) number _r, and
takes bits [_s,_s+31] of the result.*/
#define QR_FIXMUL(_a,_b,_r,_s) ((int)(((_a)*(long long)(_b)+(_r))>>(_s)))
/*Multiplies 32-bit numbers _a and _b, adds (possibly 64-bit) number _r, and
gives all 64 bits of the result.*/
#define QR_EXTMUL(_a,_b,_r) ((_a)*(long long)(_b)+(_r))
unsigned qr_isqrt(unsigned _val);
unsigned qr_ihypot(int _x,int _y);
int qr_ilog(unsigned _val);
#endif