Alien-FreeImage
view release on metacpan or search on metacpan
src/Source/FreeImageToolkit/BSplineRotate.cpp view on Meta::CPAN
// COVERED CODE IS PROVIDED UNDER THIS LICENSE ON AN "AS IS" BASIS, WITHOUT WARRANTY
// OF ANY KIND, EITHER EXPRESSED OR IMPLIED, INCLUDING, WITHOUT LIMITATION, WARRANTIES
// THAT THE COVERED CODE IS FREE OF DEFECTS, MERCHANTABLE, FIT FOR A PARTICULAR PURPOSE
// OR NON-INFRINGING. THE ENTIRE RISK AS TO THE QUALITY AND PERFORMANCE OF THE COVERED
// CODE IS WITH YOU. SHOULD ANY COVERED CODE PROVE DEFECTIVE IN ANY RESPECT, YOU (NOT
// THE INITIAL DEVELOPER OR ANY OTHER CONTRIBUTOR) ASSUME THE COST OF ANY NECESSARY
// SERVICING, REPAIR OR CORRECTION. THIS DISCLAIMER OF WARRANTY CONSTITUTES AN ESSENTIAL
// PART OF THIS LICENSE. NO USE OF ANY COVERED CODE IS AUTHORIZED HEREUNDER EXCEPT UNDER
// THIS DISCLAIMER.
//
// Use at your own risk!
// ==========================================================
/*
==========================================================
This code was taken and adapted from the following reference :
[1] Philippe Thévenaz, Spline interpolation, a C source code
implementation. http://bigwww.epfl.ch/thevenaz/
It implements ideas described in the following papers :
[2] Unser M., Splines: A Perfect Fit for Signal and Image Processing.
IEEE Signal Processing Magazine, vol. 16, no. 6, pp. 22-38, November 1999.
[3] Unser M., Aldroubi A., Eden M., B-Spline Signal Processing: Part I--Theory.
IEEE Transactions on Signal Processing, vol. 41, no. 2, pp. 821-832, February 1993.
[4] Unser M., Aldroubi A., Eden M., B-Spline Signal Processing: Part II--Efficient Design and Applications.
IEEE Transactions on Signal Processing, vol. 41, no. 2, pp. 834-848, February 1993.
==========================================================
*/
#include <float.h>
#include "FreeImage.h"
#include "Utilities.h"
#define PI ((double)3.14159265358979323846264338327950288419716939937510)
#define ROTATE_QUADRATIC 2L // Use B-splines of degree 2 (quadratic interpolation)
#define ROTATE_CUBIC 3L // Use B-splines of degree 3 (cubic interpolation)
#define ROTATE_QUARTIC 4L // Use B-splines of degree 4 (quartic interpolation)
#define ROTATE_QUINTIC 5L // Use B-splines of degree 5 (quintic interpolation)
/////////////////////////////////////////////////////////////////////////////////////////////////////////////
// Prototypes definition
static void ConvertToInterpolationCoefficients(double *c, long DataLength, double *z, long NbPoles, double Tolerance);
static double InitialCausalCoefficient(double *c, long DataLength, double z, double Tolerance);
static void GetColumn(double *Image, long Width, long x, double *Line, long Height);
static void GetRow(double *Image, long y, double *Line, long Width);
static double InitialAntiCausalCoefficient(double *c, long DataLength, double z);
static void PutColumn(double *Image, long Width, long x, double *Line, long Height);
static void PutRow(double *Image, long y, double *Line, long Width);
static bool SamplesToCoefficients(double *Image, long Width, long Height, long spline_degree);
static double InterpolatedValue(double *Bcoeff, long Width, long Height, double x, double y, long spline_degree);
static FIBITMAP * Rotate8Bit(FIBITMAP *dib, double angle, double x_shift, double y_shift, double x_origin, double y_origin, long spline_degree, BOOL use_mask);
/////////////////////////////////////////////////////////////////////////////////////////////////////////////
// Coefficients routines
/**
ConvertToInterpolationCoefficients
@param c Input samples --> output coefficients
@param DataLength Number of samples or coefficients
@param z Poles
@param NbPoles Number of poles
@param Tolerance Admissible relative error
*/
static void
ConvertToInterpolationCoefficients(double *c, long DataLength, double *z, long NbPoles, double Tolerance) {
double Lambda = 1;
long n, k;
// special case required by mirror boundaries
if(DataLength == 1L) {
return;
}
// compute the overall gain
for(k = 0L; k < NbPoles; k++) {
Lambda = Lambda * (1.0 - z[k]) * (1.0 - 1.0 / z[k]);
}
// apply the gain
for (n = 0L; n < DataLength; n++) {
c[n] *= Lambda;
}
// loop over all poles
for (k = 0L; k < NbPoles; k++) {
// causal initialization
c[0] = InitialCausalCoefficient(c, DataLength, z[k], Tolerance);
// causal recursion
for (n = 1L; n < DataLength; n++) {
c[n] += z[k] * c[n - 1L];
}
// anticausal initialization
c[DataLength - 1L] = InitialAntiCausalCoefficient(c, DataLength, z[k]);
// anticausal recursion
for (n = DataLength - 2L; 0 <= n; n--) {
c[n] = z[k] * (c[n + 1L] - c[n]);
}
}
}
/**
InitialCausalCoefficient
@param c Coefficients
@param DataLength Number of coefficients
@param z Actual pole
@param Tolerance Admissible relative error
@return
*/
static double
InitialCausalCoefficient(double *c, long DataLength, double z, double Tolerance) {
double Sum, zn, z2n, iz;
long n, Horizon;
src/Source/FreeImageToolkit/BSplineRotate.cpp view on Meta::CPAN
w2 -= w;
w4 = w2 * w2;
w -= 1.0 / 2.0;
t = w2 * (w2 - 3.0);
yWeight[0] = (1.0 / 24.0) * (1.0 / 5.0 + w2 + w4) - yWeight[5];
t0 = (1.0 / 24.0) * (w2 * (w2 - 5.0) + 46.0 / 5.0);
t1 = (-1.0 / 12.0) * w * (t + 4.0);
yWeight[2] = t0 + t1;
yWeight[3] = t0 - t1;
t0 = (1.0 / 16.0) * (9.0 / 5.0 - t);
t1 = (1.0 / 24.0) * w * (w4 - w2 - 5.0);
yWeight[1] = t0 + t1;
yWeight[4] = t0 - t1;
break;
default:
// Invalid spline degree
return 0;
}
// apply the mirror boundary conditions
for(k = 0; k <= spline_degree; k++) {
xIndex[k] = (Width == 1L) ? (0L) : ((xIndex[k] < 0L) ?
(-xIndex[k] - Width2 * ((-xIndex[k]) / Width2))
: (xIndex[k] - Width2 * (xIndex[k] / Width2)));
if (Width <= xIndex[k]) {
xIndex[k] = Width2 - xIndex[k];
}
yIndex[k] = (Height == 1L) ? (0L) : ((yIndex[k] < 0L) ?
(-yIndex[k] - Height2 * ((-yIndex[k]) / Height2))
: (yIndex[k] - Height2 * (yIndex[k] / Height2)));
if (Height <= yIndex[k]) {
yIndex[k] = Height2 - yIndex[k];
}
}
// perform interpolation
interpolated = 0.0;
for(j = 0; j <= spline_degree; j++) {
p = Bcoeff + (yIndex[j] * Width);
w = 0.0;
for(i = 0; i <= spline_degree; i++) {
w += xWeight[i] * p[xIndex[i]];
}
interpolated += yWeight[j] * w;
}
return interpolated;
}
/////////////////////////////////////////////////////////////////////////////////////////////////////////////
// FreeImage implementation
/**
Image translation and rotation using B-Splines.
@param dib Input 8-bit greyscale image
@param angle Output image rotation in degree
@param x_shift Output image horizontal shift
@param y_shift Output image vertical shift
@param x_origin Output origin of the x-axis
@param y_origin Output origin of the y-axis
@param spline_degree Output degree of the B-spline model
@param use_mask Whether or not to mask the image
@return Returns the translated & rotated dib if successful, returns NULL otherwise
*/
static FIBITMAP *
Rotate8Bit(FIBITMAP *dib, double angle, double x_shift, double y_shift, double x_origin, double y_origin, long spline_degree, BOOL use_mask) {
double *ImageRasterArray;
double p;
double a11, a12, a21, a22;
double x0, y0, x1, y1;
long x, y;
long spline;
bool bResult;
int bpp = FreeImage_GetBPP(dib);
if(bpp != 8) {
return NULL;
}
int width = FreeImage_GetWidth(dib);
int height = FreeImage_GetHeight(dib);
switch(spline_degree) {
case ROTATE_QUADRATIC:
spline = 2L; // Use splines of degree 2 (quadratic interpolation)
break;
case ROTATE_CUBIC:
spline = 3L; // Use splines of degree 3 (cubic interpolation)
break;
case ROTATE_QUARTIC:
spline = 4L; // Use splines of degree 4 (quartic interpolation)
break;
case ROTATE_QUINTIC:
spline = 5L; // Use splines of degree 5 (quintic interpolation)
break;
default:
spline = 3L;
}
// allocate output image
FIBITMAP *dst = FreeImage_Allocate(width, height, bpp);
if(!dst)
return NULL;
// buid a grey scale palette
RGBQUAD *pal = FreeImage_GetPalette(dst);
for(int i = 0; i < 256; i++) {
pal[i].rgbRed = pal[i].rgbGreen = pal[i].rgbBlue = (BYTE)i;
}
// allocate a temporary array
ImageRasterArray = (double*)malloc(width * height * sizeof(double));
if(!ImageRasterArray) {
FreeImage_Unload(dst);
return NULL;
}
// copy data samples
for(y = 0; y < height; y++) {
double *pImage = &ImageRasterArray[y*width];
BYTE *src_bits = FreeImage_GetScanLine(dib, height-1-y);
for(x = 0; x < width; x++) {
pImage[x] = (double)src_bits[x];
}
}
// convert between a representation based on image samples
// and a representation based on image B-spline coefficients
bResult = SamplesToCoefficients(ImageRasterArray, width, height, spline);
if(!bResult) {
FreeImage_Unload(dst);
free(ImageRasterArray);
return NULL;
}
// prepare the geometry
angle *= PI / 180.0;
a11 = cos(angle);
a12 = -sin(angle);
a21 = sin(angle);
a22 = cos(angle);
x0 = a11 * (x_shift + x_origin) + a12 * (y_shift + y_origin);
y0 = a21 * (x_shift + x_origin) + a22 * (y_shift + y_origin);
x_shift = x_origin - x0;
y_shift = y_origin - y0;
// visit all pixels of the output image and assign their value
for(y = 0; y < height; y++) {
BYTE *dst_bits = FreeImage_GetScanLine(dst, height-1-y);
x0 = a12 * (double)y + x_shift;
y0 = a22 * (double)y + y_shift;
for(x = 0; x < width; x++) {
x1 = x0 + a11 * (double)x;
y1 = y0 + a21 * (double)x;
if(use_mask) {
if((x1 <= -0.5) || (((double)width - 0.5) <= x1) || (y1 <= -0.5) || (((double)height - 0.5) <= y1)) {
p = 0;
}
else {
p = (double)InterpolatedValue(ImageRasterArray, width, height, x1, y1, spline);
}
}
else {
p = (double)InterpolatedValue(ImageRasterArray, width, height, x1, y1, spline);
}
// clamp and convert to BYTE
dst_bits[x] = (BYTE)MIN(MAX((int)0, (int)(p + 0.5)), (int)255);
}
}
// free working array and return
free(ImageRasterArray);
return dst;
}
/**
Image rotation using a 3rd order (cubic) B-Splines.
@param dib Input dib (8, 24 or 32-bit)
@param angle Output image rotation
@param x_shift Output image horizontal shift
@param y_shift Output image vertical shift
@param x_origin Output origin of the x-axis
@param y_origin Output origin of the y-axis
@param use_mask Whether or not to mask the image
@return Returns the translated & rotated dib if successful, returns NULL otherwise
*/
FIBITMAP * DLL_CALLCONV
FreeImage_RotateEx(FIBITMAP *dib, double angle, double x_shift, double y_shift, double x_origin, double y_origin, BOOL use_mask) {
int x, y, bpp;
int channel, nb_channels;
BYTE *src_bits, *dst_bits;
FIBITMAP *src8 = NULL, *dst8 = NULL, *dst = NULL;
if(!FreeImage_HasPixels(dib)) return NULL;
try {
bpp = FreeImage_GetBPP(dib);
if(bpp == 8) {
FIBITMAP *dst_8 = Rotate8Bit(dib, angle, x_shift, y_shift, x_origin, y_origin, ROTATE_CUBIC, use_mask);
if(dst_8) {
// copy metadata from src to dst
FreeImage_CloneMetadata(dst_8, dib);
}
return dst_8;
}
if((bpp == 24) || (bpp == 32)) {
// allocate dst image
int width = FreeImage_GetWidth(dib);
int height = FreeImage_GetHeight(dib);
if( bpp == 24 ) {
dst = FreeImage_Allocate(width, height, bpp, FI_RGBA_RED_MASK, FI_RGBA_GREEN_MASK, FI_RGBA_BLUE_MASK);
} else {
dst = FreeImage_Allocate(width, height, bpp, FI_RGBA_RED_MASK, FI_RGBA_GREEN_MASK, FI_RGBA_BLUE_MASK);
}
if(!dst) throw(1);
// allocate a temporary 8-bit dib (no need to build a palette)
src8 = FreeImage_Allocate(width, height, 8);
if(!src8) throw(1);
// process each channel separately
// -------------------------------
nb_channels = (bpp / 8);
for(channel = 0; channel < nb_channels; channel++) {
// extract channel from source dib
for(y = 0; y < height; y++) {
src_bits = FreeImage_GetScanLine(dib, y);
dst_bits = FreeImage_GetScanLine(src8, y);
for(x = 0; x < width; x++) {
dst_bits[x] = src_bits[channel];
src_bits += nb_channels;
}
}
// process channel
dst8 = Rotate8Bit(src8, angle, x_shift, y_shift, x_origin, y_origin, ROTATE_CUBIC, use_mask);
if(!dst8) throw(1);
// insert channel to destination dib
for(y = 0; y < height; y++) {
src_bits = FreeImage_GetScanLine(dst8, y);
dst_bits = FreeImage_GetScanLine(dst, y);
for(x = 0; x < width; x++) {
dst_bits[channel] = src_bits[x];
dst_bits += nb_channels;
}
}
FreeImage_Unload(dst8);
}
FreeImage_Unload(src8);
// copy metadata from src to dst
FreeImage_CloneMetadata(dst, dib);
return dst;
}
} catch(int) {
if(src8) FreeImage_Unload(src8);
if(dst8) FreeImage_Unload(dst8);
if(dst) FreeImage_Unload(dst);
}
return NULL;
}
( run in 0.568 second using v1.01-cache-2.11-cpan-5623c5533a1 )