}    /* Cases where nbr. left Lambdas = nbr. right Lambdas */    soi = 0;    while(stop2--) {        vL = vect + soi;   /* first Lambda to be read */        j = nTilde;        do {            *(vL) += (*(vG)*(*(lcPtr++)));   /* lift Lambda (Lambda + lifting value) */            vL    += stepIncr;               /* jump to next Lambda coefficient */        } while(--j);   /* use all nTilde lifting coefficients */        /* Go to next Gamma coefficient */        vG  += stepIncr;        /* Move start point for the Lambdas */        soi += stepIncr;    }    /* Cases where nbr. left Lambdas = nbr. right Lambdas */    vect += (soi - stepIncr);   /* first Lambda is always in this place */    while(stop3--) {        vL = vect;   /* position Lambda pointer */        j = nTilde;        do {            *(vL) += (*(vG)*(*(lcPtr++)));   /* lift Lambda (Lambda + lifting value) */            vL    += stepIncr;               /* jump to next Lambda coefficient */        } while(--j);   /* use all nTilde lifting coefficients */        /* Go to next Gamma coefficient */        vG += stepIncr;    }}/* CODE FOR THE IN-PLACE FLWT (AVERAGE INTERPOLATION CASE, N ODD) *//* * FLWT2D_Haar function: calculates the Haar wavelet coefficients for *                       a two dimensional signal. This case corresponds *                       to average interpolation with N=1 and wavelets *                       with nTilde = 1. */voidFLWT2D_Haar( Matrix Data, const int width, const int height,                          const int levels, const boolean inverse ){    int i, k, x, y, /* counters */        n,          /* maximum number of levels */        nX, nY,     /* number of iterations in X and Y directions */        step,       /* step size of current level */        max_step;   /* maximum and initial step for the direct and */                    /* inverse transforms, respectively */    Flt maxN;       /* maximum number of vanishing moments */    /* Number of iterations */    nX = (width==1) ? 0 :        (int)ceil((Flt)logBaseN((Flt)(width),2.0));    nX = (nX<0) ? 0 : MIN(nX,levels);    nY = (height==1) ? 0 :        (int)ceil((Flt)logBaseN((Flt)(height),2.0));    nY = (nX<0) ? 0 : MIN(nY,levels);    n = MAX(nX,nY);    if ( !inverse ) {    /* Forward transform */        /* Calculate wavelet coefficients */        max_step = (n==0) ? 0 : 1<<n;   /* maximum step */        for ( step = 1 ; step < max_step ; step <<= 1 ) {            if (nX-- > 0) {                /* Apply forward transform to the rows */                for ( y=0 ; y<height ; y+=step ) {                    FLWT_Haar ( &Data[y][0], (long)width, 1, (long)step );                }            }                if (nY-- > 0) {                /* Apply forward transform to the columns */                for ( x=0 ; x<width ; x+=step ) {                    FLWT_Haar ( &Data[0][x], (long)height, (long)width,                                  (long)step );                }            }        }    } else {   /* Inverse transform */        /* Find original dataset */        max_step = (n==0) ? 0 : 1<<(n-1);   /* initial step */        for ( step = max_step ; step >= 1 ; step >>= 1 ) {            if (n <= nY) {                /* Apply inverse transform to the columns */                for ( x=0 ; x<width ; x+=step ) {                    IFLWT_Haar ( &Data[0][x], (long)height, (long)width,                                   (long)step );                }            }                if (n <= nX) {                /* Apply inverse transform to the rows */                for ( y=0 ; y<height ; y+=step ) {                    IFLWT_Haar ( &Data[y][0], (long)width, 1, (long)step );                }            }                n--;   /* one less level to go */        }    }}/* * FLWT1D_Haar function: calculates the Haar wavelet coefficients for *                       a one dimensional signal. This case corresponds *                       to average interpolation with N=1 and wavelets *                       with nTilde = 1. */voidFLWT1D_Haar ( Vector Data, const int width,                           const int levels, const boolean inverse ){    int i, k, x, s, /* counters */        nX;         /* number of iterations in X and Y directions */    Flt maxN;       /* maximum number of vanishing moments */    /* Number of iterations */    nX = (width==1) ? 0 :        (int)ceil((Flt)logBaseN((Flt)(width),2.0));    nX = (nX<0) ? 0 : MIN(nX,levels);    if ( !inverse ) {    /* Forward transform */        /* Calculate wavelet coefficients */        for ( s = 0 ; s < nX ; s++ ) {            /* Apply forward transform to the vector */            FLWT_Haar ( Data, (long)width, 1, (long)(1<<s) );        }    } else {   /* Inverse transform */        /* Find original dataset */        for ( s = nX-1 ; s >=0 ; s-- ) {            /* Apply inverse transform to the vector */            IFLWT_Haar ( Data, (long)width, 1, (long)(1<<s) );         }    }}/* * FLWT_Haar function: finds the Haar and Lambda coefficients for *                     the given 1D signal. Boundary cases are not *                     considered yet. Thus, all signals must be *                     of even lengths. */voidFLWT_Haar ( Vector vect, const long size, const long incr,                         const long step ){    register Vector lambdaPtr,   /* pointer to Lambda coefficients */                    gammaPtr;    /* pointer to Gamma coefficients */    register long nL,            /* number of Lambda coefficients */                  nG,            /* number of Gamma coefficients */                  stepIncr;      /* step size between coefficients of the same type */    long len;                    /* number of coefficients at current level */    Flt  lAve,                   /* average of Lambdas at current level */         lSum = 0;               /* cumulative of Lambdas before last one */    /************************************************/    /* Calculate values of some important variables */    /************************************************/    len      = CEIL(size, step);   /* number of coefficients at current level */    stepIncr = (step*incr) << 1;   /* step size betweeen coefficients */    /**********************************/    /* Calculate number of iterations */    /**********************************/    nG = (len >> 1);    nL = len - nG;    /****************************/    /* Calculate Lambda average */    /****************************/    if ( ODD(len) ) {        long n = len, s = stepIncr >> 1;        Flt  t = 0;        Flt *l = vect;        while (n--) {            t += *l;            l += s;        }        lAve = t/(Flt)len;   /* average we've to maintain */    }    /**********************************************************/    /* Calculate Gamma (Haar wavelet) coefficients (odd guys) */    /* and coarser Lambda coefficients (even guys)            */    /**********************************************************/    gammaPtr = vect + (stepIncr >> 1);   /* second coefficient is Gamma */    lambdaPtr = vect;                    /* first Lambda to be read */    while (nG--) {        /* New Gamma (Haar) coefficient */        *(gammaPtr) -= *(lambdaPtr);	/* New Lambda coefficient */        *(lambdaPtr) += (*(gammaPtr)/(Flt)2);        lSum += *(lambdaPtr);        /* Go to next Gamma coefficient */        gammaPtr += stepIncr;        /* Go to next Lambda coefficient */        lambdaPtr += stepIncr;    }    /*****************************************************/    /* If ODD(len), then we have one additional case for */    /* the Lambda calculations. This is a boundary case  */    /* and, therefore, has different filter values.      */    /*****************************************************/    if ( ODD(len) ) {        *(lambdaPtr) = lAve*(Flt)nL - lSum;        lSum += *(lambdaPtr);    }}/* * IFLWT_Haar function: applies the inverse transform to find Lambda  *                      coefficients for the given 1D signal. Boundary  *                      cases are not considered yet. Thus, all signals  *                      must be of even lengths. */voidIFLWT_Haar( Vector vect, const long size, const long incr,                         const long step ){    register Vector lambdaPtr,   /* pointer to Lambda coefficients */                    gammaPtr;    /* pointer to Gamma coefficients */    register long nL,            /* number of Lambda coefficients */                  nG,            /* number of Gamma coefficients */                  stepIncr;      /* step size between coefficients of the same type */    long len;                    /* number of coefficients at current level */    Flt  lAve,                   /* average of Lambdas at current level */         lSum = 0;               /* cumulative of Lambdas before last one */    /************************************************/    /* Calculate values of some important variables */    /************************************************/    len      = CEIL(size, step);   /* number of coefficients at current level */    stepIncr = (step*incr) << 1;   /* step size betweeen coefficients */    /**********************************/    /* Calculate number of iterations */    /**********************************/    nG = (len >> 1);    nL = len - nG;    /****************************/    /* Calculate Lambda average */    /****************************/    if ( ODD(len) ) {        long n = nL, s = stepIncr;        Flt  t = 0;        Flt *l = vect;        while (n--) {            t += *l;            l += s;        }        lAve = t/(Flt)nL;   /* average we've to maintain */    }    /**********************************************************/    /* Calculate Gamma (Haar wavelet) coefficients (odd guys) */    /* and coarser Lambda coefficients (even guys)            */    /**********************************************************/    gammaPtr = vect + (stepIncr >> 1);   /* second coefficient is Gamma */    lambdaPtr = vect;                    /* first Lambda to be read */    while (nG--) {	/* New Lambda coefficient */        *(lambdaPtr) -= (*(gammaPtr)/(Flt)2);        lSum += *(lambdaPtr);	/* New Gamma coefficient */	*(gammaPtr) += *(lambdaPtr);        lSum += *(gammaPtr);        /* Go to next Gamma coefficient */        gammaPtr += stepIncr;        /* Go to next Lambda coefficient */        lambdaPtr += stepIncr;    }        /*****************************************************/    /* If ODD(len), then we have one additional case for */    /* the Lambda calculations. This is a boundary case  */    /* and, therefore, has different filter values.      */    /*****************************************************/    if ( ODD(len) ) {        *(lambdaPtr) = lAve*(Flt)len - lSum;