/*************************************************************************
* zcross.c - zero-crossings
*
* $Id: zcross.c,v 1.1 2005/12/02 22:22:14 ziad Exp $
*
* CopyrightŠINRIA 2000
*
* DESCRIPTION:
*
*
* AUTHOR:
* Gregoire Malandain (greg@sophia.inria.fr)
*
* CREATION DATE:
* Tue Nov 28 10:00:36 MET 2000
*
* Copyright Gregoire Malandain, INRIA
*
* ADDITIONS, CHANGES
*
*/
#include <zcross.h>
static int _VERBOSE_ = 0;
#define EXIT_ON_FAILURE 0
#define EXIT_ON_SUCCESS 1
void ZeroCrossings_verbose ( )
{
_VERBOSE_ = 1;
Recbuffer_verbose ( );
}
void ZeroCrossings_noverbose ( )
{
_VERBOSE_ = 0;
Recbuffer_noverbose ( );
}
/*
* `sign' of the zero-crossings
* if > 0, zero-crossings are chosen in the `positive' region
* ie are points M such that I(M) > 0 et I(M+v) <= 0
* if < 0 zero-crossings are chosen in the `negative' region
* ie are points M such that I(M) < 0 et I(M+v) >= 0
*/
static int sign_ZeroCrossing = 1;
void ZeroCrossings_Are_Positive()
{
sign_ZeroCrossing = 1;
}
void ZeroCrossings_Are_Negative()
{
sign_ZeroCrossing = -1;
}
int Extract_ZeroCrossing_2D ( void *bufferIn, bufferType typeIn,
void *bufferOut, bufferType typeOut, int *bufferDims )
{
if ( sign_ZeroCrossing > 0 )
return( Extract_PositiveZeroCrossing_2D( bufferIn, typeIn, bufferOut,
typeOut, bufferDims ) );
return( Extract_NegativeZeroCrossing_2D( bufferIn, typeIn, bufferOut,
typeOut, bufferDims ) );
}
/*************************************************************
ZERO-CROSSINGS
*/
int Extract_PositiveZeroCrossing_2D ( void *bufferIn,
bufferType typeIn,
void *bufferOut,
bufferType typeOut,
int *bufferDims )
{
/* les passage par zeros sont definis par
I(M) > 0 et I(M+v) <= 0
on marque M
Dans un contexte de detection de contours (ex laplacien)
cela marque dans les zones sombres (par rapport a clair)
*/
char *proc="Extract_PositiveZeroCrossing_2D";
int x, y, z, iz, iy;
int dx = bufferDims[0];
int dxy = bufferDims[0] * bufferDims[1];
/*
* We check the buffers' dimensions.
*/
if ( (bufferDims[0] <= 0) || (bufferDims[1] <= 0) || (bufferDims[2] <= 0) ) {
if ( _VERBOSE_ > 0 )
fprintf( stderr, " Fatal error in %s: improper buffer's dimension.\n", proc );
return( EXIT_ON_FAILURE );
}
if ( bufferIn == bufferOut ) {
if ( _VERBOSE_ > 0 )
fprintf( stderr, " Error in %s: input buffer should not be equal to output.\n", proc );
return( EXIT_ON_FAILURE );
}
/*
* // the (bufferDims[1]-1) first rows
* for ( y=0, iy=0; y<bufferDims[1]-1; y++, iy+=dx )
*
* // the (dx-1) first points of the row
* for ( x=0; x<dx-1; x++ )
* // we mark (x,y) if (x,y) > 0 && ((x+1,y) <= 0 || (x,y+1) <= 0)
* // we mark (x+1,y) if (x,y) <= 0 && (x+1,y) > 0
* // we mark (x,y+1) if (x,y) <= 0 && (x,y+1) > 0
*
* // the last point of the row
* // we mark (x,y) if (x,y) > 0 && (x,y+1) <= 0
* // we mark (x,y+1) if (x,y) <= 0 && (x,y+1) > 0
*
* // the last row
* for ( x=0; x<dx-1; x++ )
* // we mark (x,y) if (x,y) > 0 && (x+1,y) <= 0
* // we mark (x+1,y) if (x,y) <= 0 && (x+1,y) > 0
*
*/
#define _POSITIVE_ZERO_CROSSINGS_( TYPE ) { \
TYPE *resBuf = (TYPE*)bufferOut; \
iz = bufferDims[2]*bufferDims[1]*bufferDims[0]; \
for ( x=0; x<iz; x++ ) resBuf[x] = 0; \
for ( z=0, iz=0; z<bufferDims[2]; z++, iz+=dxy ) { \
for ( y=0, iy=0; y<bufferDims[1]-1; y++, iy+=dx ) { \
for ( x=0; x<dx-1; x++ ) { \
if ( theBuf[iz+iy+x] > 0 ) { \
if ( theBuf[iz+iy+x+1] <= 0 || theBuf[iz+iy+x+dx] <= 0 ) resBuf[iz+iy+x] = 1; \
} else { \
if ( theBuf[iz+iy+x+1] > 0 ) resBuf[iz+iy+x+1] = 1; \
if ( theBuf[iz+iy+x+dx] > 0 ) resBuf[iz+iy+x+dx] = 1; \
} \
} \
if ( theBuf[iz+iy+x] > 0 ) { \
if ( theBuf[iz+iy+x+dx] <= 0 ) resBuf[iz+iy+x] = 1; \
} else { \
if ( theBuf[iz+iy+x+dx] > 0 ) resBuf[iz+iy+x+dx] = 1; \
} \
} \
for ( x=0; x<dx-1; x++ ) { \
if ( theBuf[iz+iy+x] > 0 ) { \
if ( theBuf[iz+iy+x+1] <= 0 ) resBuf[iz+iy+x] = 1; \
} else { \
if ( theBuf[iz+iy+x+1] > 0 ) resBuf[iz+iy+x+1] = 1; \
} \
} \
} \
}
switch( typeIn ) {
case FLOAT :
{
r32 *theBuf = (r32*)bufferIn;
switch( typeOut ) {
case UCHAR :
_POSITIVE_ZERO_CROSSINGS_( u8 )
break;
case FLOAT :
_POSITIVE_ZERO_CROSSINGS_( r32 )
break;
default :
if ( _VERBOSE_ > 0 )
fprintf( stderr, " Error in %s: such output type not handled.\n", proc );
return( EXIT_ON_FAILURE );
}
} /* typeIn == FLOAT */
break;
default :
if ( _VERBOSE_ > 0 )
fprintf( stderr, " Error in %s: such input type not handled.\n", proc );
return( EXIT_ON_FAILURE );
}
return( EXIT_ON_SUCCESS );
}
int Extract_NegativeZeroCrossing_2D ( void *bufferIn,
bufferType typeIn,
void *bufferOut,
bufferType typeOut,
int *bufferDims )
{
/* les passage par zeros sont definis par
I(M) < 0 et I(M+v) >= 0
on marque M
Dans un contexte de detection de contours (ex laplacien)
cela marque dans les zones sombres (par rapport a clair)
*/
char *proc="Extract_NegativeZeroCrossing_2D";
int x, y, z, iz, iy;
int dx = bufferDims[0];
int dxy = bufferDims[0] * bufferDims[1];
/*
* We check the buffers' dimensions.
*/
if ( (bufferDims[0] <= 0) || (bufferDims[1] <= 0) || (bufferDims[2] <= 0) ) {
if ( _VERBOSE_ > 0 )
fprintf( stderr, " Fatal error in %s: improper buffer's dimension.\n", proc );
return( EXIT_ON_FAILURE );
}
if ( bufferIn == bufferOut ) {
if ( _VERBOSE_ > 0 )
fprintf( stderr, " Error in %s: input buffer should not be equal to output.\n", proc );
return( EXIT_ON_FAILURE );
}
/*
* // the (bufferDims[1]-1) first rows
* for ( y=0, iy=0; y<bufferDims[1]-1; y++, iy+=dx )
*
* // the (dx-1) first points of the row
* for ( x=0; x<dx-1; x++ )
* // we mark (x,y) if (x,y) < 0 && ((x+1,y) >= 0 || (x,y+1) >= 0)
* // we mark (x+1,y) if (x,y) >= 0 && (x+1,y) < 0
* // we mark (x,y+1) if (x,y) >= 0 && (x,y+1) < 0
*
* // the last point of the row
* // we mark (x,y) if (x,y) < 0 && (x,y+1) >= 0
* // we mark (x,y+1) if (x,y) >= 0 && (x,y+1) < 0
*
* // the last row
* for ( x=0; x<dx-1; x++ )
* // we mark (x,y) if (x,y) < 0 && (x+1,y) >= 0
* // we mark (x+1,y) if (x,y) >= 0 && (x+1,y) < 0
*
*/
#define _NEGATIVE_ZERO_CROSSINGS_( TYPE ) { \
TYPE *resBuf = (TYPE*)bufferOut; \
iz = bufferDims[2]*bufferDims[1]*bufferDims[0]; \
for ( x=0; x<iz; x++ ) resBuf[x] = 0; \
for ( z=0, iz=0; z<bufferDims[2]; z++, iz+=dxy ) { \
for ( y=0, iy=0; y<bufferDims[1]-1; y++, iy+=dx ) { \
for ( x=0; x<dx-1; x++ ) { \
if ( theBuf[iz+iy+x] < 0 ) { \
if ( theBuf[iz+iy+x+1] >= 0 || theBuf[iz+iy+x+dx] >= 0 ) resBuf[iz+iy+x] = 1; \
} else { \
if ( theBuf[iz+iy+x+1] < 0 ) resBuf[iz+iy+x+1] = 1; \
if ( theBuf[iz+iy+x+dx] < 0 ) resBuf[iz+iy+x+dx] = 1; \
} \
} \
if ( theBuf[iz+iy+x] < 0 ) { \
if ( theBuf[iz+iy+x+dx] >= 0 ) resBuf[iz+iy+x] = 1; \
} else { \
if ( theBuf[iz+iy+x+dx] < 0 ) resBuf[iz+iy+x+dx] = 1; \
} \
} \
for ( x=0; x<dx-1; x++ ) { \
if ( theBuf[iz+iy+x] < 0 ) { \
if ( theBuf[iz+iy+x+1] >= 0 ) resBuf[iz+iy+x] = 1; \
} else { \
if ( theBuf[iz+iy+x+1] < 0 ) resBuf[iz+iy+x+1] = 1; \
} \
} \
} \
}
switch( typeIn ) {
case FLOAT :
{
r32 *theBuf = (r32*)bufferIn;
switch( typeOut ) {
case UCHAR :
_NEGATIVE_ZERO_CROSSINGS_( u8 )
break;
case FLOAT :
_NEGATIVE_ZERO_CROSSINGS_( r32 )
break;
default :
if ( _VERBOSE_ > 0 )
fprintf( stderr, " Error in %s: such output type not handled.\n", proc );
return( EXIT_ON_FAILURE );
}
}
break;
default :
if ( _VERBOSE_ > 0 )
fprintf( stderr, " Error in %s: such input type not handled.\n", proc );
return( EXIT_ON_FAILURE );
}
return( EXIT_ON_SUCCESS );
}
int Mask_With_Image( void *bufferIn, bufferType typeIn,
void *bufferMask, bufferType typeMask,
void *bufferOut, bufferType typeOut,
int *bufferDims )
{
char *proc = "Mask_With_Image";
int i, v;
/*
* We check the buffers' dimensions.
*/
if ( (bufferDims[0] <= 0) || (bufferDims[1] <= 0) || (bufferDims[2] <= 0) ) {
if ( _VERBOSE_ > 0 )
fprintf( stderr, " Fatal error in %s: improper buffer's dimension.\n", proc );
return( EXIT_ON_FAILURE );
}
v = bufferDims[0] * bufferDims[1] * bufferDims[2];
if ( typeIn != typeOut ) {
if ( _VERBOSE_ > 0 )
fprintf( stderr, " Fatal error in %s: buffers in and out should have the same type.\n", proc );
return( EXIT_ON_FAILURE );
}
#define _MASK_( type ) { \
type *theBuf = (type*)bufferIn; \
type *resBuf = (type*)bufferOut; \
for (i=0; i<v; i++ ) { \
resBuf[i] = ( theMask[i] > 0 ) ? theBuf[i] : 0.0; \
} \
}
switch( typeMask ) {
case UCHAR :
{
u8 *theMask = (u8*)bufferMask;
switch( typeIn ) {
case FLOAT :
_MASK_( r32 )
break;
case DOUBLE :
_MASK_( r64 )
break;
default :
if ( _VERBOSE_ > 0 )
fprintf( stderr, " Error in %s: such output type not handled.\n", proc );
return( EXIT_ON_FAILURE );
}
} /* mask: end of UCHAR */
break;
case FLOAT :
{
r32 *theMask = (r32*)bufferMask;
switch( typeIn ) {
case FLOAT :
_MASK_( r32 )
break;
case DOUBLE :
_MASK_( r64 )
break;
default :
if ( _VERBOSE_ > 0 )
fprintf( stderr, " Error in %s: such output type not handled.\n", proc );
return( EXIT_ON_FAILURE );
}
} /* mask: end of FLOAT */
break;
default :
if ( _VERBOSE_ > 0 )
fprintf( stderr, " Error in %s: such mask type not handled.\n", proc );
return( EXIT_ON_FAILURE );
}
return( EXIT_ON_SUCCESS );
}
/*
*
*
* Edge detection
* - laplacian
* - hessian
*
*
*/
static void GradientModulus2D( float *gradient_modulus,
float *derivative_along_X,
float *derivative_along_Y,
int length )
{
register int i;
register float *norme = gradient_modulus;
register float *gx = derivative_along_X;
register float *gy = derivative_along_Y;
for ( i=0; i<length; i++, norme++, gx++, gy++ )
*norme = sqrt( (*gx)*(*gx) + (*gy)*(*gy) );
}
static void GradientModulus3D( float *gradient_modulus,
float *derivative_along_X,
float *derivative_along_Y,
float *derivative_along_Z,
int length )
{
register int i;
register float *norme = gradient_modulus;
register float *gx = derivative_along_X;
register float *gy = derivative_along_Y;
register float *gz = derivative_along_Z;
for ( i=0; i<length; i++, norme++, gx++, gy++, gz++ )
*norme = sqrt( (*gx)*(*gx) + (*gy)*(*gy) + (*gz)*(*gz) );
}
int Gradient_On_Laplacian_ZeroCrossings_2D ( void *bufferIn, bufferType typeIn,
void *bufferOut, bufferType typeOut,
int *bufferDims, int *borderLengths,
float *filterCoefs,
recursiveFilterType filterType )
{
char *proc = "Gradient_On_Laplacian_ZeroCrossings_2D";
float *theGR = NULL;
float *theXX = NULL;
float *theYY = NULL;
derivativeOrder Xderiv[3] = { DERIVATIVE_1_EDGES, SMOOTHING, NODERIVATIVE };
derivativeOrder Yderiv[3] = { SMOOTHING, DERIVATIVE_1_EDGES, NODERIVATIVE };
derivativeOrder XXderiv[3] = { DERIVATIVE_2, SMOOTHING, NODERIVATIVE };
derivativeOrder YYderiv[3] = { SMOOTHING, DERIVATIVE_2, NODERIVATIVE };
int sliceDims[3];
int z, i, dimxXdimy;
void *sliceOut = NULL;
/*
* We check the buffers' dimensions.
*/
if ( (bufferDims[0] <= 0) || (bufferDims[1] <= 0) || (bufferDims[2] <= 0) ) {
if ( _VERBOSE_ > 0 )
fprintf( stderr, " Fatal error in %s: improper buffer's dimension.\n", proc );
return( EXIT_ON_FAILURE );
}
/*
* test of the coefficients
*/
if ( (filterCoefs[0] < 0.0) || (filterCoefs[1] < 0.0) ||
(filterCoefs[2] < 0.0) ) {
if ( _VERBOSE_ > 0 )
fprintf( stderr, " Error in %s: negative coefficient's value.\n", proc );
return( EXIT_ON_FAILURE );
}
/*
*
*/
dimxXdimy = bufferDims[0] * bufferDims[1];
sliceDims[0] = bufferDims[0];
sliceDims[1] = bufferDims[1];
sliceDims[2] = 1;
if ( typeOut == FLOAT ) {
theXX = (float*)malloc( 2 * dimxXdimy * sizeof( float ) );
} else {
theXX = (float*)malloc( 3 * dimxXdimy * sizeof( float ) );
}
if ( theXX == NULL ) {
if ( _VERBOSE_ > 0 ) {
fprintf( stderr, " Fatal error in %s:", proc );
fprintf( stderr, " unable to allocate auxiliary buffer.\n" );
}
return( EXIT_ON_FAILURE );
}
theGR = theXX;
theGR += dimxXdimy;
if ( typeOut != FLOAT ) {
theYY = theGR;
theYY += dimxXdimy;
}
for ( z=0; z<bufferDims[2]; z++ ) {
if ( typeOut == FLOAT ) {
theYY = ((float*)bufferOut) + z * dimxXdimy;
}
if ( RecursiveFilterOnBuffer( bufferIn, typeIn, theXX, FLOAT,
sliceDims, borderLengths,
Xderiv, filterCoefs, filterType ) == 0 ) {
if ( _VERBOSE_ > 0 ) {
fprintf( stderr, " Fatal error in %s:", proc );
fprintf( stderr, " unable to compute X derivative.\n" );
}
free( theXX );
return( EXIT_ON_FAILURE );
}
if ( RecursiveFilterOnBuffer( bufferIn, typeIn, theGR, FLOAT,
sliceDims, borderLengths,
Yderiv, filterCoefs, filterType ) == 0 ) {
if ( _VERBOSE_ > 0 ) {
fprintf( stderr, " Fatal error in %s:", proc );
fprintf( stderr, " unable to compute Y derivative.\n" );
}
free( theXX );
return( EXIT_ON_FAILURE );
}
GradientModulus2D( theGR, theGR, theXX, dimxXdimy );
if ( RecursiveFilterOnBuffer( bufferIn, typeIn, theXX, FLOAT,
sliceDims, borderLengths,
XXderiv, filterCoefs, filterType ) == 0 ) {
if ( _VERBOSE_ > 0 ) {
fprintf( stderr, " Fatal error in %s:", proc );
fprintf( stderr, " unable to compute X^2 derivative.\n" );
}
free( theXX );
return( EXIT_ON_FAILURE );
}
if ( RecursiveFilterOnBuffer( bufferIn, typeIn, theYY, FLOAT,
sliceDims, borderLengths,
YYderiv, filterCoefs, filterType ) == 0 ) {
if ( _VERBOSE_ > 0 ) {
fprintf( stderr, " Fatal error in %s:", proc );
fprintf( stderr, " unable to compute Y^2 derivative.\n" );
}
free( theXX );
return( EXIT_ON_FAILURE );
}
/*
* theYY = laplacian
*/
for ( i=0; i<dimxXdimy; i++ ) theYY[i] += theXX[i];
/*
* theXX = zero-crossings
*/
if ( Extract_ZeroCrossing_2D( theYY, FLOAT, theXX, FLOAT, sliceDims ) == 0 ) {
if ( _VERBOSE_ > 0 ) {
fprintf( stderr, " Fatal error in %s:", proc );
fprintf( stderr, " unable to compute zero crossing.\n" );
}
free( theXX );
return( EXIT_ON_FAILURE );
}
if ( Mask_With_Image( theGR, FLOAT, theXX, FLOAT, theYY, FLOAT, sliceDims ) == 0 ) {
if ( _VERBOSE_ > 0 ) {
fprintf( stderr, " Fatal error in %s:", proc );
fprintf( stderr, " unable to mask with zero crossing.\n" );
}
free( theXX );
return( EXIT_ON_FAILURE );
}
if ( typeOut != FLOAT ) {
switch ( typeOut ) {
case UCHAR :
sliceOut = (((u8*)bufferOut) + z * dimxXdimy);
break;
case SCHAR :
sliceOut = (((s8*)bufferOut) + z * dimxXdimy);
break;
case SSHORT :
sliceOut = (((s16*)bufferOut) + z * dimxXdimy);
break;
case DOUBLE :
sliceOut = (((r64*)bufferOut) + z * dimxXdimy);
break;
default :
if ( _VERBOSE_ > 0 )
fprintf( stderr, " Error in %s: such output type not handled.\n", proc );
free( theXX );
return( EXIT_ON_FAILURE );
}
ConvertBuffer( theYY, FLOAT, sliceOut, typeOut, dimxXdimy );
}
}
return( EXIT_ON_SUCCESS );
}
int Gradient_On_GradientHessianGradient_ZeroCrossings_2D ( void *bufferIn,
bufferType typeIn,
void *bufferOut,
bufferType typeOut,
int *bufferDims,
int *borderLengths,
float *filterCoefs,
recursiveFilterType filterType )
{
char *proc = "Gradient_On_GradientHessianGradient_ZeroCrossings_2D";
float *theXX = NULL;
float *theYY = NULL;
float *theXY = NULL;
float *theX = NULL;
float *theY = NULL;
derivativeOrder Xsmooth[3] = { SMOOTHING, NODERIVATIVE, NODERIVATIVE };
derivativeOrder Yderiv[3] = { NODERIVATIVE, DERIVATIVE_1_EDGES, NODERIVATIVE };
derivativeOrder YYderiv[3] = { NODERIVATIVE, DERIVATIVE_2, NODERIVATIVE };
derivativeOrder Ysmooth[3] = { NODERIVATIVE, SMOOTHING, NODERIVATIVE };
derivativeOrder Xderiv[3] = { DERIVATIVE_1_EDGES, NODERIVATIVE, NODERIVATIVE };
derivativeOrder XXderiv[3] = { DERIVATIVE_2, NODERIVATIVE, NODERIVATIVE };
derivativeOrder XYderiv[3] = { DERIVATIVE_1, DERIVATIVE_1, NODERIVATIVE };
int sliceDims[3];
int z, i, dimxXdimy;
void *sliceOut = NULL;
double gx, gy;
/*
* We check the buffers' dimensions.
*/
if ( (bufferDims[0] <= 0) || (bufferDims[1] <= 0) || (bufferDims[2] <= 0) ) {
if ( _VERBOSE_ > 0 )
fprintf( stderr, " Fatal error in %s: improper buffer's dimension.\n", proc );
return( EXIT_ON_FAILURE );
}
/*
* test of the coefficients
*/
if ( (filterCoefs[0] < 0.0) || (filterCoefs[1] < 0.0) ||
(filterCoefs[2] < 0.0) ) {
if ( _VERBOSE_ > 0 )
fprintf( stderr, " Error in %s: negative coefficient's value.\n", proc );
return( EXIT_ON_FAILURE );
}
/*
*
*/
dimxXdimy = bufferDims[0] * bufferDims[1];
sliceDims[0] = bufferDims[0];
sliceDims[1] = bufferDims[1];
sliceDims[2] = 1;
if ( typeOut == FLOAT ) {
theXX = (float*)malloc( 4 * dimxXdimy * sizeof( float ) );
} else {
theXX = (float*)malloc( 5 * dimxXdimy * sizeof( float ) );
}
if ( theXX == NULL ) {
if ( _VERBOSE_ > 0 ) {
fprintf( stderr, " Fatal error in %s:", proc );
fprintf( stderr, " unable to allocate auxiliary buffer.\n" );
}
return( EXIT_ON_FAILURE );
}
theX = theY = theYY = theXX;
theYY += dimxXdimy;
theX += 2*dimxXdimy;
theY += 3*dimxXdimy;
if ( typeOut != FLOAT ) {
theXY = theXX;
theXY += 4*dimxXdimy;
}
for ( z=0; z<bufferDims[2]; z++ ) {
if ( typeOut == FLOAT ) {
theXY = ((float*)bufferOut) + z * dimxXdimy;
}
if ( RecursiveFilterOnBuffer( bufferIn, typeIn, theX, FLOAT,
sliceDims, borderLengths,
Ysmooth, filterCoefs, filterType ) == 0 ) {
if ( _VERBOSE_ > 0 ) {
fprintf( stderr, " Fatal error in %s:", proc );
fprintf( stderr, " unable to compute Y^0 derivative.\n" );
}
free( theXX );
return( EXIT_ON_FAILURE );
}
if ( RecursiveFilterOnBuffer( bufferIn, typeIn, theY, FLOAT,
sliceDims, borderLengths,
Xsmooth, filterCoefs, filterType ) == 0 ) {
if ( _VERBOSE_ > 0 ) {
fprintf( stderr, " Fatal error in %s:", proc );
fprintf( stderr, " unable to compute X^0 derivative.\n" );
}
free( theXX );
return( EXIT_ON_FAILURE );
}
if ( RecursiveFilterOnBuffer( bufferIn, typeIn, theXY, FLOAT,
sliceDims, borderLengths,
XYderiv, filterCoefs, filterType ) == 0 ) {
if ( _VERBOSE_ > 0 ) {
fprintf( stderr, " Fatal error in %s:", proc );
fprintf( stderr, " unable to compute X^1Y^1 derivative.\n" );
}
free( theXX );
return( EXIT_ON_FAILURE );
}
if ( RecursiveFilterOnBuffer( theX, FLOAT, theXX, FLOAT,
sliceDims, borderLengths,
XXderiv, filterCoefs, filterType ) == 0 ) {
if ( _VERBOSE_ > 0 ) {
fprintf( stderr, " Fatal error in %s:", proc );
fprintf( stderr, " unable to compute X^2 derivative.\n" );
}
free( theXX );
return( EXIT_ON_FAILURE );
}
if ( RecursiveFilterOnBuffer( theY, FLOAT, theYY, FLOAT,
sliceDims, borderLengths,
YYderiv, filterCoefs, filterType ) == 0 ) {
if ( _VERBOSE_ > 0 ) {
fprintf( stderr, " Fatal error in %s:", proc );
fprintf( stderr, " unable to compute Y^2 derivative.\n" );
}
free( theXX );
return( EXIT_ON_FAILURE );
}
if ( RecursiveFilterOnBuffer( theX, FLOAT, theX, FLOAT,
sliceDims, borderLengths,
Xderiv, filterCoefs, filterType ) == 0 ) {
if ( _VERBOSE_ > 0 ) {
fprintf( stderr, " Fatal error in %s:", proc );
fprintf( stderr, " unable to compute X^1 derivative.\n" );
}
free( theXX );
return( EXIT_ON_FAILURE );
}
if ( RecursiveFilterOnBuffer( theY, FLOAT, theY, FLOAT,
sliceDims, borderLengths,
Yderiv, filterCoefs, filterType ) == 0 ) {
if ( _VERBOSE_ > 0 ) {
fprintf( stderr, " Fatal error in %s:", proc );
fprintf( stderr, " unable to compute Y^1 derivative.\n" );
}
free( theXX );
return( EXIT_ON_FAILURE );
}
/*
* theXY = gradient . Hessian * gradient
*/
for ( i=0; i<dimxXdimy; i++ ) {
gx = theX[i];
gy = theY[i];
theXY[i] = gx * ( theXX[i] * gx + theXY[i] * gy ) +
gy * ( theXY[i] * gx + theYY[i] * gy );
}
/*
* theYY = gradient modulus
*/
GradientModulus2D( theYY, theX, theY, dimxXdimy );
/*
* theXX = zero-crossings
*/
if ( Extract_ZeroCrossing_2D( theXY, FLOAT, theXX, FLOAT, sliceDims ) == 0 ) {
if ( _VERBOSE_ > 0 ) {
fprintf( stderr, " Fatal error in %s:", proc );
fprintf( stderr, " unable to compute zero crossing.\n" );
}
free( theXX );
return( EXIT_ON_FAILURE );
}
if ( Mask_With_Image( theYY, FLOAT, theXX, FLOAT, theXY, FLOAT, sliceDims ) == 0 ) {
if ( _VERBOSE_ > 0 ) {
fprintf( stderr, " Fatal error in %s:", proc );
fprintf( stderr, " unable to mask with zero crossing.\n" );
}
free( theXX );
return( EXIT_ON_FAILURE );
}
if ( typeOut != FLOAT ) {
switch ( typeOut ) {
case UCHAR :
sliceOut = (((u8*)bufferOut) + z * dimxXdimy);
break;
case SCHAR :
sliceOut = (((s8*)bufferOut) + z * dimxXdimy);
break;
case SSHORT :
sliceOut = (((s16*)bufferOut) + z * dimxXdimy);
break;
case DOUBLE :
sliceOut = (((r64*)bufferOut) + z * dimxXdimy);
break;
default :
if ( _VERBOSE_ > 0 )
fprintf( stderr, " Error in %s: such output type not handled.\n", proc );
free( theXX );
return( EXIT_ON_FAILURE );
}
ConvertBuffer( theXY, FLOAT, sliceOut, typeOut, dimxXdimy );
}
}
return( EXIT_ON_SUCCESS );
}
syntax highlighted by Code2HTML, v. 0.9.1