suffixTree.c

Go to the documentation of this file.

/* Copyright 2013 Jorge Merlino

   This file is part of Context.

   Context is free software: you can redistribute it and/or modify
   it under the terms of the GNU General Public License as published by
   the Free Software Foundation, either version 3 of the License, or
   (at your option) any later version.

   Context is distributed in the hope that it will be useful,
   but WITHOUT ANY WARRANTY; without even the implied warranty of
   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
   GNU General Public License for more details.

   You should have received a copy of the GNU General Public License
   along with Context.  If not, see <http://www.gnu.org/licenses/>.
*/
 
#include <stdlib.h>
#include "suffixTree.h"
#include "stack.h"
#include "alpha.h"
#include "text.h"
#include "gammaFunc.h"
#include "spacedef.h"
#include "debug.h"

#define INFINITY -3

#define BOTTOM -2

#define ROOT -1

#define GET_RIGHT(R, T)\
        if (T->child) {\
          R = T->child->left-1;\
        }\
        else {\
          R = INFINITY;\
        }\


static void canonize (suffixTree_t s, Uint k, Uint p, suffixTree_t *retS, Uint *retK) {
  Uint kPrime, pPrime;
  suffixTree_t sPrime;

  if (p == -1 || p < k) {
    *retS = s;
    *retK = k;
  }
  else {
    sPrime = s;
    if (s->left == BOTTOM) {
      sPrime = s->child; /* ROOT is the t-transition of BOTTOM for every t */
    }
    else {
      for (sPrime=sPrime->child; sPrime && text[sPrime->left]!=text[k]; sPrime=sPrime->sibling);
    }
    kPrime = sPrime->left;
    GET_RIGHT(pPrime, sPrime);
    while ((k <= p) && (pPrime - kPrime <= p - k)) {
      k = k + pPrime - kPrime +1;
      s = sPrime;
      if (k <= p) {
    for (sPrime=sPrime->child; sPrime && text[sPrime->left]!=text[k]; sPrime=sPrime->sibling);
    kPrime = sPrime->left;
    GET_RIGHT(pPrime, sPrime);
      }
    }
    *retS = s;
    *retK = k;
  }
}


static BOOL testAndSplit(suffixTree_t s, Uint k, Uint p, suffixTree_t *r) {
  suffixTree_t sPrime = s, new;
  Uint kPrime;

  if (p != -1 && k <= p) {
    for (sPrime=sPrime->child; sPrime && text[sPrime->left]!=text[k]; sPrime=sPrime->sibling);
    kPrime = sPrime->left;
    if ((p+1 < textlen) && (text[kPrime+p-k+1] == text[p+1])) {
      *r = s;
      return True;
    }
    else {
      /* split transition */
      CALLOC(new, struct suffixTree, 1);
      new->left = sPrime->left;
      new->child = sPrime;
      new->sibling = sPrime->sibling;
      new->parent = sPrime->parent;

      sPrime->left = kPrime + p - k + 1;
      sPrime->parent = new;
      sPrime->sibling = NULL;

      /* fix the parent */
      if (new->parent->child == sPrime) {
    new->parent->child = new;
      }
      else {
    for (s=new->parent->child; s->sibling != sPrime; s=s->sibling);
    s->sibling = new;
      }

      *r = new;
      return False;
    }
  }
  else { /* s is explicit */
    *r = s;
    if (s->left == BOTTOM) {
      return True;
    }
    if (p+1 < textlen) {
      for (sPrime=sPrime->child; sPrime && text[sPrime->left]!=text[p+1]; sPrime=sPrime->sibling);
    }
    else {
      sPrime=NULL;
    }
    return sPrime != NULL;
  }
}


static void addChild (suffixTree_t s, Uint i) {
  suffixTree_t new;

  CALLOC(new, struct suffixTree, 1);
  new->left = i;
  new->parent = s;
  
  if (!s->child) {
    s->child = new;
  }
  else {
    for (s=s->child; s->sibling; s=s->sibling);
    s->sibling = new;
  }
}


static void update (suffixTree_t s, Uint k, Uint i, suffixTree_t *retS, Uint *retK) {
  suffixTree_t oldr=NULL, r;

  while (!testAndSplit(s, k, i-1, &r)) {
    /* add new leaf */
    addChild(r, i);
    if (oldr) {
      /* create suffix link */
      oldr->suffix = r;
    }
    oldr = r;
    canonize(s->suffix, k, i-1, &s, &k);
  }

  if (oldr && oldr->left!=ROOT) { 
    oldr->suffix = r;
  }
  *retS = s;
  *retK = k;
}


static void pruneSubTree(suffixTree_t tree) {
  Uint stacktop=0, stackalloc=0, *stack = NULL, treePtr;

  PUSHNODE((Uint)tree);
  while(NOTSTACKEMPTY) {
    POPNODE(treePtr);
    tree = (suffixTree_t)treePtr;
    if (tree->sibling) {
      PUSHNODE((Uint)tree->sibling);
    }
    if (tree->child) {
      PUSHNODE((Uint)tree->child);
    }
    freeSuffixTree(tree);
  }
  FREE(stack);
}

#ifdef DEBUG 

static void printNode(suffixTree_t tree) {
  Uint right;

  if (tree->left == BOTTOM) {
    /* nothing */
  }
  else if (tree->left == ROOT) {
    printf("root %p", tree);
  }
  else {
    GET_RIGHT(right, tree);
    if (right == INFINITY) right = textlen;
    printf("%ld - %ld - %p", tree->left, right, tree);
  }
}


static void printRec(suffixTree_t tree, Uint level) {
  int i;

  if (level > 0) {
    for (i=1; i<level; i++) {
      printf("-");
    }
    printNode (tree);
    printf("\n");
  }
  
  level++;
  if (tree->child) {
    for (tree=tree->child; tree; tree=tree->sibling) {
      printRec(tree, level);
    }
  }
}


static int treeSize(suffixTree_t tree) {
  int sum = 1;

  if (tree->child) {
    for (tree=tree->child; tree; tree=tree->sibling) {
      sum += treeSize(tree);
    }
  }
  return sum;
}

#endif


suffixTree_t initSuffixTree() {
  suffixTree_t tree;

  CALLOC(tree, struct suffixTree, 1);
  CALLOC(tree->child, struct suffixTree, 1);

  tree->left = BOTTOM;
  
  tree->child->left = ROOT;
  tree->child->suffix = tree;
  tree->child->parent = tree;

  return tree;
}


void freeSuffixTree(suffixTree_t tree) {
  FREE(tree);
}


void buildSuffixTree(suffixTree_t tree) {
  Uint k=0, i=0;
  suffixTree_t s = tree->child; /* ROOT */

  while (i <= textlen) {
    update(s, k, i, &s, &k);
    canonize(s, k, i, &s, &k);
    i++;
    if (i%1000 == 0)  printf("%ld/%ld\r", i, textlen);
  }
  printf("%ld/%ld\n", textlen, textlen);
  DEBUGCODE(printf("Original tree size: %d\n", treeSize(tree)));
}


void pruneSuffixTree(suffixTree_t tree) {
  Uint stacktop=0, stackalloc=0, *stack = NULL, treePtr, length, i, right;
  suffixTree_t child;
  double est, auxx;
  Uint * distinct;

  tree = tree->child; /* ROOT */
  PUSHNODE((Uint)tree);
  PUSHNODE(0);
  while(NOTSTACKEMPTY) {
    POPNODE(length);
    POPNODE(treePtr);
    tree = (suffixTree_t)treePtr;
    if (!tree->child) { /* is a leaf */
      tree->stats = allocStatistics();
      if (tree->left > length) { /* this is not the full string prefix */
    i = GETINDEX(tree->left-length-1);
    tree->stats->count[i] = 1;
    tree->stats->symbols[tree->stats->symbolCount++] = i;
      }
      /*tree->stats->cost = log2Alpha();*/
    }
    else { /* is a branching node */
      for (child=tree->child; child && child->stats; child=child->sibling);
      if (child) { /* there is one child to evaluate */
    PUSHNODE((Uint)tree);
    PUSHNODE(length);
    PUSHNODE((Uint)child);
    if (tree->left != ROOT) {
      GET_RIGHT(right, tree);
      PUSHNODE(length + right - tree->left + 1);
    }
    else {
      PUSHNODE(0);
    }
      } 
      else { /* all children are already evaluated */
    tree->stats = allocStatistics();
    CALLOC(distinct, Uint, alphasize);

    for (child=tree->child; child; child=child->sibling) {
      for (i=0; i<child->stats->symbolCount; i++) {
        if (tree->stats->count[child->stats->symbols[i]] == 0) {
          tree->stats->symbols[tree->stats->symbolCount++] = child->stats->symbols[i];
        } 
        tree->stats->count[child->stats->symbols[i]] += child->stats->count[child->stats->symbols[i]];
        distinct[child->stats->symbols[i]]++;
      }
      tree->stats->cost += child->stats->cost;
      freeStatistics(child->stats);
    }

    GET_RIGHT(right, tree);
    tree->stats->cost += hAlpha() * alphasize * (right - tree->left + 1); 

    auxx = escapeCost(tree->stats, distinct);
    tree->stats->cost += auxx;

    free(distinct);

    est = nodeCost(tree->stats);

    if (est <= tree->stats->cost) { /* we have to prune */
      tree->stats->cost = est;
      pruneSubTree(tree->child);
      tree->child = NULL;
    }
      }
    }
  }
  freeStatistics(tree->stats);
  FREE(stack);
}


fsmTree_t fsmSuffixTree(suffixTree_t tree) {
  Uint stacktop=0, stackalloc=0, *stack = NULL, sfxPtr, fsmPtr, pos, right;
  fsmTree_t ret = initFsmTree(), fsmNode;

  if (tree->child->child) { /* if root has children */
    PUSHNODE((Uint)tree->child->child); 
    PUSHNODE((Uint)ret); /* ROOT */
    while(NOTSTACKEMPTY) {
      POPNODE(fsmPtr);
      POPNODE(sfxPtr);
      fsmNode = (fsmTree_t)fsmPtr;
      tree = (suffixTree_t)sfxPtr;
      if (tree->left != textlen) { /* ignore $ leaves */
    pos = GETINDEX(tree->left);   
    fsmNode->children[pos] = initFsmTree();
    fsmNode->children[pos]->parent = fsmNode; 
    fsmNode->children[pos]->left = tree->left; 
    GET_RIGHT(right, tree);
    if (right == INFINITY) {
      fsmNode->children[pos]->right = tree->left; /* shorten leaf */
    }
    else {
      fsmNode->children[pos]->right = right; 
    }
    if (fsmNode->left != ROOT) {
      fsmNode->children[pos]->length = fsmNode->length + fsmNode->right - fsmNode->left + 1;
    }

    if (tree->sibling) {
      PUSHNODE((Uint)tree->sibling);
      PUSHNODE(fsmPtr);
    }
    if (tree->child) {
      PUSHNODE((Uint)tree->child);
      PUSHNODE((Uint)fsmNode->children[pos]);
    }
      }
      FREE(tree);
    }
    FREE(stack);
  }
  return ret;
}

#ifdef DEBUG

void printSuffixTree(const suffixTree_t tree) {
  printRec(tree, 0);
}

#endif