 /*
  * Khoros: $Id: lvfht.c,v 1.1 1991/05/10 15:41:54 khoros Exp $
  */

#if !defined(lint) && !defined(SABER)
static char rcsid[] = "Khoros: $Id: lvfht.c,v 1.1 1991/05/10 15:41:54 khoros Exp $";
#endif

 /*
  * $Log: lvfht.c,v $
 * Revision 1.1  1991/05/10  15:41:54  khoros
 * Initial revision
 *
  */ 

/*
 *----------------------------------------------------------------------
 *
 * Copyright 1991, University of New Mexico.  All rights reserved.
 * Permission to copy and modify this software and its documen-
 * tation only for internal use in your organization is hereby
 * granted, provided that this notice is retained thereon and
 * on all copies.  UNM makes no representations as too the sui-
 * tability and operability of this software for any purpose.
 * It is provided "as is" without express or implied warranty.
 * 
 * UNM DISCLAIMS ALL WARRANTIES WITH REGARD TO THIS SOFTWARE,
 * INCLUDING ALL IMPLIED WARRANTIES OF MERCHANTABILITY AND FIT-
 * NESS.  IN NO EVENT SHALL UNM BE LIABLE FOR ANY SPECIAL,
 * INDIRECT OR CONSEQUENTIAL DAMAGES OR ANY OTHER DAMAGES WHAT-
 * SOEVER RESULTING FROM LOSS OF USE, DATA OR PROFITS, WHETHER
 * IN AN ACTION OF CONTRACT, NEGLIGENCE OR OTHER TORTIOUS
 * ACTION, ARISING OUT OF OR IN CONNECTION WITH THE USE OR PER-
 * FORMANCE OF THIS SOFTWARE.
 * 
 * No other rights, including for example, the right to redis-
 * tribute this software and its documentation or the right to
 * prepare derivative works, are granted unless specifically
 * provided in a separate license agreement.
 *---------------------------------------------------------------------
 */

#include "unmcopyright.h"        /* Copyright 1991 by UNM */

/*>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>  <<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<
 >>>>
 >>>>         File Name: lvfht.c
 >>>>
 >>>>      Program Name: vfht
 >>>>
 >>>> Date Last Updated: Tue Mar  5 22:41:28 1991 
 >>>>
 >>>>          Routines: lvfht - the library call for vfht
 >>>>
 >>>>>>>>>>>>>>>>>>>>>>>>>>>>>>   <<<<<<<<<<<<<<<<<<<<<<<<<<<<<<*/


#include "vinclude.h"


/* -library_includes */
#define CHECKLIBINPUT(program, img) { \
        if (!propertype (program, img, VFF_TYP_FLOAT, FALSE)) {               \
            fprintf (stderr, "error in lvfht: input image must ");          \
            fprintf (stderr, "be VFF_TYP_FLOAT: returning with an error\n");  \
            return (0);                                                       \
        }                                                                     \
        if (!proper_loc_type (program, img, VFF_LOC_IMPLICIT, FALSE)) {       \
            fprintf (stderr, "error in lvfht: input image must ");          \
            fprintf (stderr,"be VFF_LOC_IMPLICIT: returning with an error\n");\
            return (0);                                                       \
        }                                                                     \
        if (!proper_map_type (program, img, VFF_MAPTYP_NONE, FALSE)) {        \
            fprintf (stderr, "error in lvfht: input image must ");          \
            fprintf (stderr,"be VFF_MAPTYP_NONE: returning with an error\n"); \
            return (0);                                                       \
        }                                                                     \
        if (!proper_map_enable (program, img, VFF_MAP_OPTIONAL, FALSE)) {     \
            fprintf (stderr, "error in lvfht: input image must ");          \
            fprintf (stderr,"be VFF_MAP_OPTIONAL: returning with an error\n");\
            return (0);                                                       \
        }                                                                     \
        if (!proper_map_scheme (program, img, VFF_MS_NONE, FALSE)) {          \
            fprintf (stderr, "error in lvfht: input image must ");          \
            fprintf (stderr,"be VFF_MS_NONE: returning with an error\n");     \
            return (0);                                                       \
        }                                                                     \
        if (!proper_num_images (program, img, 1, FALSE)) {                    \
            fprintf (stderr, "error in lvfht: input image must ");          \
            fprintf(stderr,"contain only 1 image: returning with an error\n");\
            return (0);                                                       \
        }                                                                     \
        if (!proper_num_bands (program, img, 1, FALSE)) {                     \
            fprintf (stderr, "error in lvfht: input image must ");          \
            fprintf(stderr,"contain only 1 band: returning with an error\n"); \
            return (0);                                                       \
        } }

static void fht1d();
/* -library_includes_end */


/****************************************************************
*
* Routine Name: lvfht - library call for vfht
*
* Purpose:
*    
*    Take the 2-D Fast Hadamard Transform of an image.
*    
*    
* Input:
*    
*         1.  img - a pointer to an xvimage structure  containing  the
*         image to be transformed
*    
*         2.  inv_flag - if set this flag indicates  that  an  inverse
*         transform should be performed, otherwise a forward transform
*         is executed.
*    
*    
* Output:
*    
*         1.  img - the input image is overwritten and the  output  is
*         returned  in  its place.  Be careful not to overwrite impor-
*         tant data.
*    
*         lvfht returns 1 upon success and 0 on failure.
*    
*    
*
* Written By: Per Lysne
*    
*    
****************************************************************/


/* -library_def */
int lvfht (img, inv_flag)
struct xvimage *img;
int inv_flag;
/* -library_def_end */

/* -library_code */
{
    int N, M, i, j;
    float *data, *column;
    char *program="lvfht";

    /* Check for a bad input pointer. */
    if (img==NULL) {
        fprintf (stderr, "error in lvfht: bad input xvimage pointer: ");
        fprintf (stderr, "returning false\n");
        return (0);
    }

    /* Check for a correct type of input image. */
    CHECKLIBINPUT (program, img);

    /* Check for a legal size of input image. */
    if (!powtwo(img->row_size)) {
        fprintf (stderr,"error in lvfht: row length of the input image must");
        fprintf (stderr, " be a power of two: exiting with an error\n");
        return (0);
    }
    if (!powtwo(img->col_size)) {
        fprintf (stderr,"error in lvfht: col length of the input image must");
        fprintf (stderr, " be a power of two: exiting with an error\n");
        return (0);
    }

    /* N is the number of samples. */
    M = img->row_size;
    N = img->col_size;

    column = (float *) calloc(N, sizeof(float));
    if (column==NULL) {
        fprintf (stderr, "error in lvfht: could not allocate column: ");
        fprintf (stderr, "exiting false\n");
        return (0);
    }

    /* Use data so I don't have to say img->imagedata all the time. */
    data=(float *) img->imagedata;

    /* Do the FHT across the rows. */
    for (i=0; i<N; i++)
        fht1d(&data[i*M], M, inv_flag);

    /* Do the FWT down the columns. */
    for (j=0; j<M; j++ ) {
        for (i=0; i<N; i++)
            column[i]=data[i*M+j];
        fht1d(column, N, inv_flag);
        for (i=0; i<N; i++)
            data[i*M+j]=column[i];
    }

    return(1);
}



static void fht1d (data, N, inv_flag)
float *data;
int N;
int inv_flag;
{
    int i, ip, j, k, l, le;
    float temp;

    /*
     * This loop reorders the input data.
     */
    j=0;
    for (i=0; i<N-1; i++) {
        if (i<j) {
            temp=data[j];
            data[j]=data[i];
            data[i]=temp;
        }
        k=N/2;
        while (k<=j) {
            j=j-k;
            k=k/2;
        }
        j=j+k;
    }
    
    /*
     * This loop does the Hadamard Transform.
     */
    l=1;
    while (l<N) {
        le=l;
        l=l*2;
        for (j=0; j<le; j++) {
            for (i=j; i<N; i=i+l) {
                ip=i+le;
                temp=data[ip];
                data[ip]=data[i]-temp;
                data[i]=data[i]+temp;
            }
        }
    }

    /* Scale the output of a forward transform. */
    if (inv_flag)
        for (i=0; i<N; i++)
            data[i]=data[i]/N;

    return;
}
/* -library_code_end */
