"Fossies" - the Fresh Open Source Software Archive  

Source code changes of the file "Libtmp/Image2D/resample.c" between
PDL-2.080.tar.gz and PDL-2.081.tar.gz

About: PDL (Perl Data Language) aims to turn perl into an efficient numerical language for scientific computing (similar to IDL and MatLab).

[ To the main PDL source changes report ]

resample.c  (PDL-2.080):resample.c  (PDL-2.081)
skipping to change at line 30 skipping to change at line 30
@return x to the power p. @return x to the power p.
@doc @doc
This is much faster than the math function due to the integer. Some This is much faster than the math function due to the integer. Some
compilers make this optimization already, some do not. compilers make this optimization already, some do not.
p can be positive, negative or null. p can be positive, negative or null.
*/ */
/*--------------------------------------------------------------------------*/ /*--------------------------------------------------------------------------*/
double ipow(double x, int p) long double ipow(long double x, int p)
{ {
double r, recip ; long double r, recip ;
/* Get rid of trivial cases */ /* Get rid of trivial cases */
switch (p) { switch (p) {
case 0: case 0:
return 1.00 ; return 1.00 ;
case 1: case 1:
return x ; return x ;
case 2: case 2:
skipping to change at line 70 skipping to change at line 70
} }
return r; return r;
} }
/* /*
* compute the value of a 2D polynomial at a point * compute the value of a 2D polynomial at a point
* - it assumes that ncoeff is a small number, so there's * - it assumes that ncoeff is a small number, so there's
* little point in pre-calculating ipow(u,i) * little point in pre-calculating ipow(u,i)
*/ */
double long double
poly2d_compute( int ncoeff, double *c, double u, double *vpow ) { poly2d_compute( int ncoeff, long double *c, long double u, long double *vpow ) {
double out; long double out;
int i, j, k; int i, j, k;
out = 0.00; out = 0.00;
k = 0; k = 0;
for( j = 0; j < ncoeff; j++ ) { for( j = 0; j < ncoeff; j++ ) {
for( i = 0; i < ncoeff; i++ ) { for( i = 0; i < ncoeff; i++ ) {
out += c[k] * ipow( u, i ) * vpow[j]; out += c[k] * ipow( u, i ) * vpow[j];
k++; k++;
} }
} }
skipping to change at line 122 skipping to change at line 122
@param nn Number of samples in the input kernel. @param nn Number of samples in the input kernel.
@return void @return void
@doc @doc
Bring back a hyperbolic tangent kernel from Fourier to normal space. Do Bring back a hyperbolic tangent kernel from Fourier to normal space. Do
not try to understand the implementation and DO NOT MODIFY THIS FUNCTION. not try to understand the implementation and DO NOT MODIFY THIS FUNCTION.
*/ */
/*--------------------------------------------------------------------------*/ /*--------------------------------------------------------------------------*/
#define KERNEL_SW(a,b) tempr=(a);(a)=(b);(b)=tempr #define KERNEL_SW(a,b) tempr=(a);(a)=(b);(b)=tempr
static void reverse_tanh_kernel(double * data, int nn) static void reverse_tanh_kernel(long double * data, int nn)
{ {
unsigned long n, unsigned long n,
mmax, mmax,
m, m,
i, j, i, j,
istep ; istep ;
double wtemp, long double wtemp,
wr, wr,
wpr, wpr,
wpi, wpi,
wi, wi,
theta; theta;
double tempr, long double tempr,
tempi; tempi;
n = (unsigned long)nn << 1; n = (unsigned long)nn << 1;
j = 1; j = 1;
for (i=1 ; i<n ; i+=2) { for (i=1 ; i<n ; i+=2) {
if (j > i) { if (j > i) {
KERNEL_SW(data[j-1],data[i-1]); KERNEL_SW(data[j-1],data[i-1]);
KERNEL_SW(data[j],data[i]); KERNEL_SW(data[j],data[i]);
} }
m = n >> 1; m = n >> 1;
skipping to change at line 204 skipping to change at line 204
\item It is infinitely differentiable everywhere (i.e. smooth). \item It is infinitely differentiable everywhere (i.e. smooth).
\item The transition sharpness is scalable. \item The transition sharpness is scalable.
\end{itemize} \end{itemize}
The returned array must be deallocated using free(). The returned array must be deallocated using free().
*/ */
/*--------------------------------------------------------------------------*/ /*--------------------------------------------------------------------------*/
#define hk_gen(x,s) (((tanh(s*(x+0.5))+1)/2)*((tanh(s*(-x+0.5))+1)/2)) #define hk_gen(x,s) (((tanh(s*(x+0.5))+1)/2)*((tanh(s*(-x+0.5))+1)/2))
double * generate_tanh_kernel(double steep) long double * generate_tanh_kernel(long double steep)
{ {
double * kernel ; long double * kernel ;
double * x ; long double * x ;
double width ; long double width ;
double inv_np ; long double inv_np ;
double ind ; long double ind ;
int i ; int i ;
int np ; int np ;
int samples ; int samples ;
width = (double)TABSPERPIX / 2.0 ; width = (long double)TABSPERPIX / 2.0 ;
samples = KERNEL_SAMPLES ; samples = KERNEL_SAMPLES ;
np = 32768 ; /* Hardcoded: should never be changed */ np = 32768 ; /* Hardcoded: should never be changed */
inv_np = 1.00 / (double)np ; inv_np = 1.00 / (long double)np ;
/* /*
* Generate the kernel expression in Fourier space * Generate the kernel expression in Fourier space
* with a correct frequency ordering to allow standard FT * with a correct frequency ordering to allow standard FT
*/ */
x = (double *) malloc((2*np+1)*sizeof(double)) ; x = (long double *) malloc((2*np+1)*sizeof(long double)) ;
for (i=0 ; i<np/2 ; i++) { for (i=0 ; i<np/2 ; i++) {
ind = (double)i * 2.0 * width * inv_np ; ind = (long double)i * 2.0 * width * inv_np ;
x[2*i] = hk_gen(ind, steep) ; x[2*i] = hk_gen(ind, steep) ;
x[2*i+1] = 0.00 ; x[2*i+1] = 0.00 ;
} }
for (i=np/2 ; i<np ; i++) { for (i=np/2 ; i<np ; i++) {
ind = (double)(i-np) * 2.0 * width * inv_np ; ind = (long double)(i-np) * 2.0 * width * inv_np ;
x[2*i] = hk_gen(ind, steep) ; x[2*i] = hk_gen(ind, steep) ;
x[2*i+1] = 0.00 ; x[2*i+1] = 0.00 ;
} }
/* /*
* Reverse Fourier to come back to image space * Reverse Fourier to come back to image space
*/ */
reverse_tanh_kernel(x, np) ; reverse_tanh_kernel(x, np) ;
/* /*
* Allocate and fill in returned array * Allocate and fill in returned array
*/ */
kernel = (double *) malloc(samples * sizeof(double)) ; kernel = (long double *) malloc(samples * sizeof(long double)) ;
for (i=0 ; i<samples ; i++) { for (i=0 ; i<samples ; i++) {
kernel[i] = 2.0 * width * x[2*i] * inv_np ; kernel[i] = 2.0 * width * x[2*i] * inv_np ;
} }
free(x) ; free(x) ;
return kernel ; return kernel ;
} /* generate_tanh_kernel() */ } /* generate_tanh_kernel() */
/*-------------------------------------------------------------------------*/ /*-------------------------------------------------------------------------*/
/** /**
skipping to change at line 279 skipping to change at line 279
"lanczos" & Lanczos2 kernel \\ "lanczos" & Lanczos2 kernel \\
"hamming" & Hamming kernel \\ "hamming" & Hamming kernel \\
"hann" & Hann kernel "hann" & Hann kernel
\end{tabular} \end{tabular}
The returned array of doubles is ready of use in the various re-sampling The returned array of doubles is ready of use in the various re-sampling
functions in this module. It must be deallocated using free(). functions in this module. It must be deallocated using free().
*/ */
/*--------------------------------------------------------------------------*/ /*--------------------------------------------------------------------------*/
double * long double *
generate_interpolation_kernel(char * kernel_type) generate_interpolation_kernel(char * kernel_type)
{ {
double * tab ; long double * tab ;
int i ; int i ;
double x ; long double x ;
double alpha ; long double alpha ;
double inv_norm ; long double inv_norm ;
int samples = KERNEL_SAMPLES ; int samples = KERNEL_SAMPLES ;
if (kernel_type==NULL) { if (kernel_type==NULL) {
tab = generate_interpolation_kernel("tanh") ; tab = generate_interpolation_kernel("tanh") ;
} else if (!strcmp(kernel_type, "default")) { } else if (!strcmp(kernel_type, "default")) {
tab = generate_interpolation_kernel("tanh") ; tab = generate_interpolation_kernel("tanh") ;
} else if (!strcmp(kernel_type, "sinc")) { } else if (!strcmp(kernel_type, "sinc")) {
tab = (double *) malloc(samples * sizeof(double)) ; tab = (long double *) malloc(samples * sizeof(long double)) ;
tab[0] = 1.0 ; tab[0] = 1.0 ;
tab[samples-1] = 0.0 ; tab[samples-1] = 0.0 ;
for (i=1 ; i<samples ; i++) { for (i=1 ; i<samples ; i++) {
x = (double)KERNEL_WIDTH * (double)i/(double)(samples-1) ; x = (long double)KERNEL_WIDTH * (long double)i/(long doub le)(samples-1) ;
tab[i] = sinc(x) ; tab[i] = sinc(x) ;
} }
} else if (!strcmp(kernel_type, "sinc2")) { } else if (!strcmp(kernel_type, "sinc2")) {
tab = (double *) malloc(samples * sizeof(double)) ; tab = (long double *) malloc(samples * sizeof(long double)) ;
tab[0] = 1.0 ; tab[0] = 1.0 ;
tab[samples-1] = 0.0 ; tab[samples-1] = 0.0 ;
for (i=1 ; i<samples ; i++) { for (i=1 ; i<samples ; i++) {
x = 2.0 * (double)i/(double)(samples-1) ; x = 2.0 * (long double)i/(long double)(samples-1) ;
tab[i] = sinc(x) ; tab[i] = sinc(x) ;
tab[i] *= tab[i] ; tab[i] *= tab[i] ;
} }
} else if (!strcmp(kernel_type, "lanczos")) { } else if (!strcmp(kernel_type, "lanczos")) {
tab = (double *) malloc(samples * sizeof(double)) ; tab = (long double *) malloc(samples * sizeof(long double)) ;
for (i=0 ; i<samples ; i++) { for (i=0 ; i<samples ; i++) {
x = (double)KERNEL_WIDTH * (double)i/(double)(samples-1) ; x = (long double)KERNEL_WIDTH * (long double)i/(long doub le)(samples-1) ;
if (fabs(x)<2) { if (fabs(x)<2) {
tab[i] = sinc(x) * sinc(x/2) ; tab[i] = sinc(x) * sinc(x/2) ;
} else { } else {
tab[i] = 0.00 ; tab[i] = 0.00 ;
} }
} }
} else if (!strcmp(kernel_type, "hamming")) { } else if (!strcmp(kernel_type, "hamming")) {
tab = (double *) malloc(samples * sizeof(double)) ; tab = (long double *) malloc(samples * sizeof(long double)) ;
alpha = 0.54 ; alpha = 0.54 ;
inv_norm = 1.00 / (double)(samples - 1) ; inv_norm = 1.00 / (long double)(samples - 1) ;
for (i=0 ; i<samples ; i++) { for (i=0 ; i<samples ; i++) {
x = (double)i ; x = (long double)i ;
if (i<(samples-1)/2) { if (i<(samples-1)/2) {
tab[i] = alpha + (1-alpha) * cos(2.0*PI_NUMB*x*in v_norm) ; tab[i] = alpha + (1-alpha) * cos(2.0*PI_NUMB*x*in v_norm) ;
} else { } else {
tab[i] = 0.0 ; tab[i] = 0.0 ;
} }
} }
} else if (!strcmp(kernel_type, "hann")) { } else if (!strcmp(kernel_type, "hann")) {
tab = (double *) malloc(samples * sizeof(double)) ; tab = (long double *) malloc(samples * sizeof(long double)) ;
alpha = 0.50 ; alpha = 0.50 ;
inv_norm = 1.00 / (double)(samples - 1) ; inv_norm = 1.00 / (long double)(samples - 1) ;
for (i=0 ; i<samples ; i++) { for (i=0 ; i<samples ; i++) {
x = (double)i ; x = (long double)i ;
if (i<(samples-1)/2) { if (i<(samples-1)/2) {
tab[i] = alpha + (1-alpha) * cos(2.0*PI_NUMB*x*in v_norm) ; tab[i] = alpha + (1-alpha) * cos(2.0*PI_NUMB*x*in v_norm) ;
} else { } else {
tab[i] = 0.0 ; tab[i] = 0.0 ;
} }
} }
} else if (!strcmp(kernel_type, "tanh")) { } else if (!strcmp(kernel_type, "tanh")) {
tab = generate_tanh_kernel(TANH_STEEPNESS) ; tab = generate_tanh_kernel(TANH_STEEPNESS) ;
} else { } else {
/** /**
 End of changes. 29 change blocks. 
37 lines changed or deleted 37 lines changed or added
To the C Programs section of the PDL source changes report or the top of this page
Mainly generated by pkgdiff using rfcdiff
Home  |  About  |  Features  |  All  |  Newest  |  Dox  |  Diffs  |  RSS Feeds  |  Screenshots  |  Comments  |  Imprint  |  Privacy  |  HTTP(S)