Los Alamos National Laboratory/LANL_CD.ISO/software/compres/src/make_tre.c

C/C++ Source or Header | 1992-02-11 | 6KB | 224 lines

/****************************************************************

COPYRIGHT (C) 1992 UNIVERSITY OF CALIFORNIA

***************************************************************/

#define BUCKETSIZE 16
#define NULL 0
#include <math.h>
#include "kdnode.h"

int *TreeList;

/****************************************************************/

struct KDTreeNode *MakeMedianTree(k, ncv, cbook)
int k, ncv;
float *cbook;
{
    struct KDTreeNode *TempNode;
    int i;
    float **Vecs;
    void MedianPartition();

    conv2d(cbook, &Vecs, ncv, k);

    if(!(TreeList = (int*)malloc(ncv * sizeof(int)))) {
        printf("ERROR: Not enough room for malloc\n");
        exit(-1);
    }

    for (i = 0; i < ncv; i++)
        TreeList[i] = i;

    if(!(TempNode = (struct KDTreeNode *)malloc(sizeof(struct KDTreeNode)))) {
        printf("ERROR: Not enough room for malloc\n");
        exit(-1);
    }

    MedianPartition(TempNode, TreeList, ncv, k, ncv, Vecs);
    return(TempNode);
}

/****************************************************************/

void MedianPartition(Node, ListStart, ListSize, k, ncv, Vecs)
struct KDTreeNode *Node;
int *ListStart, ListSize, k, ncv;
float **Vecs;
{
    float Median, Val, Vari, MaxVari, StDev();
    int Axis, Middle, *LList, *RList, LSize, RSize, MaxAxis;
    struct KDTreeNode *LTemp, *RTemp;
    void QuickSort();

    if (ListSize < BUCKETSIZE)  {
        Node->Left = NULL;
        Node->Right = NULL;
        Node->List = ListStart;
        Node->ListSize = ListSize;
        Node->KeyIndex = 0;
        Node->PartitionValue = 0;
        return;
    }

    else  {
        for (Axis = 0; Axis < k; Axis++)  {
            Vari = StDev (ListStart, ListSize, Axis, Vecs);
                if ((Axis == 0) || (Vari > MaxVari)) {
                    MaxVari = Vari;
                    MaxAxis = Axis;
        }
    }

    Axis = MaxAxis;

/* Find the median by sorting the list and taking the middle element */

    QuickSort(ListStart, 0, (ListSize - 1), Vecs, Axis);
    Middle = (ListSize / 2);
    Median = Vecs[Axis][ListStart[Middle]];

/* Partition the array into two lists */

/* Move the dividing point up to the last occurrence of Median */

    for ( ; Middle < ListSize; Middle++)  {
        Val = Vecs[Axis][ListStart[Middle]];
        if (Val != Median)
            break;
    }
    Middle--;

/* Determine the size of each sublist */

    LSize = Middle + 1;
    RSize = ListSize - LSize;

  /* Check to see if either sublist is of size 0: if so, make this a
   * terminal node. */

    if ((LSize == 0) || (RSize == 0))  {
        printf("DUMB CONDITION\n");
        Node->Left = NULL;
        Node->Right = NULL;
        Node->List = ListStart;
        Node->ListSize = ListSize;
        Node->KeyIndex = 0;
        Node->PartitionValue = 0;
        return;
    }

/* Get pointers to sublists */
    LList = ListStart;
    RList = &ListStart[Middle + 1];

/* Create subnodes */
    LTemp = (struct KDTreeNode*)malloc(sizeof(struct KDTreeNode));
    RTemp = (struct KDTreeNode*)malloc(sizeof(struct KDTreeNode));
    if(!LTemp || !RTemp) {
        printf("ERROR: Not enough room for malloc\n");
        exit(-1);
    }

/* Set values for this node */
    Node->Left = LTemp;
    Node->Right = RTemp;
    Node->List = NULL;
    Node->ListSize = 0;
    Node->KeyIndex = Axis;
    Node->PartitionValue = Median;

/* Partition the subnodes */
    MedianPartition (LTemp, LList, LSize, k, ncv, Vecs);
    MedianPartition (RTemp, RList, RSize, k, ncv, Vecs);

    }
}

/****************************************************************/

/* Does a Quicksort of a linear array.  The numbers in the array are index
 * numbers for the Vecs[][] array and the datum to sort by is the position
 * of each vector along axis number Axis. */

void QuickSort (Array, Start, Finish, Vecs, Axis)
int *Array, Start, Finish, Axis;
float **Vecs;
{
    int Left, Right, Temp;
    float StarterValue;

    Left = Start;
    Right = Finish;
    StarterValue = Vecs[Axis][Array[Start]];

    while (Left < Right)  {

        while (Vecs[Axis][Array[Left]] < StarterValue)
            Left++;
        while (Vecs[Axis][Array[Right]] > StarterValue)
            Right--;
        if (Left <= Right)  {
            Temp = Array[Left];
            Array[Left++] = Array[Right];
            Array[Right--] = Temp;
    }

}

    if (Start < Right)
        QuickSort (Array, Start, Right, Vecs, Axis);
    if (Left < Finish)
        QuickSort (Array, Left, Finish, Vecs, Axis);

}

/****************************************************************/

/* Accepts a list of vector indices Array and calculates the standard deviation
 * of the indicated vectors in Vecs along the axis Axis. */

float StDev (Array, ArraySize, Axis, Vecs)
int *Array, ArraySize, Axis;
float **Vecs;
{
    int Index, i;
    float Mean, Var, xx;

    Mean = 0.;
    for (i = 0; i < ArraySize; i++)  {
        Index = Array[i];
        Mean += Vecs[Axis][Index];
    }
    Mean /= (float)ArraySize;

    xx = 0.;
    for (i = 0; i < ArraySize; i++)  {
        Index = Array[i];
        Var = Vecs[Axis][Index] - Mean;
        if (Var < 0)
            Var = -Var;
        xx += Var;
    }
    xx /= (float) ArraySize;

    return xx;
}

/****************************************************************/

conv2d(cbook, Vecs, ncv, k)
int ncv, k;
float *cbook, ***Vecs;
{
    int i, j;
    float **matrix();

    *Vecs = matrix(k, ncv);
    for(i = 0; i < ncv; i++)
        for(j = 0; j < k; j++)
            (*Vecs)[j][i] = *cbook++;
}