rightext.java

A program to find frequent molecular substructures and discriminative fragments in a database of mol

/*----------------------------------------------------------------------
  File    : RightExt.java
  Contents: Molecular fragment extension management
            (depth first search/rightmost extension)
  Author  : Christian Borgelt
  History : 11.03.2002 file created as file submol.java
            11.08.2005 class Extension made abstract, second strategy
            13.08.2005 precomputation of rightmost path added
            11.04.2006 function compareBond added (for comparing bonds)
            30.04.2006 treatment of additional atoms added
            02.05.2006 functions makeCanonic and makeWord added
            10.05.2006 function makeString added, makeCanonic written
            17.05.2006 extension-specific ring function added
            18.05.2006 adapted to new field Extension.dst
            31.05.2006 function isCanonic extended, result changed
            04.06.2006 creation and check of ring extensions added
            06.06.2006 adapted to changed type of bond flags
            07.06.2006 function compareWord added
            01.07.2006 adaptation of ring extensions moved here
            10.08.2006 adapted to new class TypeMap
----------------------------------------------------------------------*/
package moss;

/*--------------------------------------------------------------------*/
public class RightExt extends Extension {
/*--------------------------------------------------------------------*/

  /* --- instance variables --- */
  private int[] path;           /* indices of bonds on rightmost path */
  private int   pbi;            /* index of next path bond */

  /*------------------------------------------------------------------*/

  public RightExt (int mode, int max, TypeMap map)
  {                             /* --- create a rightmost extension */
    super(mode, max, map);      /* call base initialization, */
    this.path = new int[256];   /* create a bond index vector, */
    this.frag = null;           /* and clear the fragment */
  }  /* RightExt() */

  /*------------------------------------------------------------------*/

  public void init (Fragment frag, Embedding base)
  {                             /* --- init. a rightmost extension */
    int  i, k, n;               /* loop variables, bond counter */
    Atom a, d;                  /* to traverse/access the atoms */
    Bond b, x;                  /* to traverse/access the bonds */
    int  v[];                   /* buffer for reallocation */

    base.markId();              /* mark the embedding in the molecule */
    this.base = base;           /* note the embedding to extend */
    this.src  = base.atoms.length -1;
    this.pbi  =  0;             /* start with the first bond of */
    this.idx  = -1;             /* the atom with the highest index */
    this.size = -1;             /* clear the ring size and flags */
    this.all  =  0;             /* (also extension type indicator) */
    if (frag == this.frag)      /* if to work on the same fragment, */
      return;                   /* no additional work is necessary */
    this.frag = frag;           /* otherwise note the new fragment */
    n = 0; b = null;            /* init. the number of bonds and */
    a = base.atoms[frag.dst];   /* get rightmost leaf (end of path) */
    while (a.mark > 0) {        /* construct rightmost path upwards */
      for (k = -1, i = a.bondcnt; --i >= 0; ) {
        x = a.bonds[i];         /* traverse the marked bonds */
        if (x.mark < 0) continue;
        d = (x.src != a) ? x.src : x.dst;
        if ((d.mark > k) && (d.mark < a.mark)) {
          k = d.mark; b = x; }  /* find the adjacent atom with the */
      }                         /* largest index smaller than own */
      if (n >= this.path.length) {
        v = new int[n +(n >> 1)];
        System.arraycopy(this.path, 0, v, 0, n);
        this.path = v;          /* if the bond index vector */
      }                         /* is too small, enlarge it */
      this.path[n++] = b.mark;  /* note the index of the next bond */
      a = base.atoms[k];        /* on the rightmost path and */
    }                           /* get the next path atom */
  }  /* init() */

  /*------------------------------------------------------------------*/

  public boolean next ()
  {                             /* --- create the next extension */
    Atom a, d, dd, p[];         /* to traverse/access the atoms */
    Bond b, bb;                 /* to traverse/access the bonds */

    /* --- continue an old extension --- */
    if (((this.mode & RING) != 0)
    &&  (this.all != 0)         /* if there is another ring flag */
    &&   this.ring())           /* and some ring is admissible, */
      return true;              /* return "extension successful" */

    /* --- find a new extension from an additional atom --- */
    b = null; d = null;         /* dummy initializations */
    p = this.base.atoms;        /* get the atoms of the base embed. */
    a = p[this.src];            /* and the current anchor atom */
    while (this.src > this.frag.dst) {  /* if at an add. atom, */
      if (++this.idx >= a.bondcnt) {    /* go to the next bond */
        a = p[--this.src];      /* if atom's last bond is processed, */
        this.idx = -1;          /* go to the preceding atom */
        continue;               /* (may be add. atom or end of path) */
      }                         /* start with its first bond */
      b = a.bonds[this.idx];    /* get the next bond of this atom */
      if (b.mark != -1)         /* if the bond is in the embedding */
        continue;               /* or excluded, it cannot be added */
      d = (a != b.src) ? b.src : b.dst;
      if ((d.mark < 0)          /* if atom is not in the embedding */
      &&  (p.length +this.frag.chns >= this.max))
        continue;               /* check whether a new atom is ok */
      if (this.store(a, b, d)) return true;
    }                           /* store the new extension */

    /* -- find a new extension from the rightmost path -- */
    while (true) {              /* find the next unprocessed bond */
      while (++this.idx >= a.bondcnt) {
        if (this.src <= 0) {    /* if atom's last bond is processed, */
          this.base.unmark();   /* check for another path bond and */
          return false;         /* if there is none (at the root), */
        }                       /* unmark the embedding and abort */
        b = this.base.bonds[this.path[this.pbi++]];
        a = b.src; d = b.dst;   /* get next path bond and its atoms */
        if (d.mark < a.mark) { a = d; d = b.src; }
        this.src = a.mark;      /* note next atom on rightmost path */
        this.idx = -1;          /* find the downward bond from it */
        while (a.bonds[++this.idx] != b);
        while (--this.idx >= 0) {  /* check for equivalent bonds */
          bb = a.bonds[this.idx];  /* among the preceding ones */
          if (bb.type != b.type) break;
          dd = (bb.src != a) ? bb.src : bb.dst;
          if (dd.type != d.type) break;
        }                       /* (equivalent processed bonds */
      }                         /*  must be considered again) */
      b = a.bonds[this.idx];    /* get the next bond of this atom */
      if (b.mark != -1)         /* if the bond is in the embedding */
        continue;               /* or excluded, it cannot be added */
      d = (a != b.src) ? b.src : b.dst;
      if ((d.mark >= 0)         /* skip bonds closing a ring that */
      && ((d.mark <  this.frag.dst)   /* do not lead to the leaf */
      ||  (a.mark >= this.frag.src))) /* or start below the source */
        continue;                     /* of the preceding bond */
      if ((d.mark < 0)          /* if atom is not in the embedding */
      &&  (p.length +this.frag.chns >= this.max))
        continue;               /* check whether a new atom is ok */
      if (this.store(a, b, d)) return true;
    }                           /* store the new extension */
  }  /* next() */

  /*------------------------------------------------------------------*/

  private boolean store (Atom a, Bond b, Atom d)
  {                             /* --- store next extension */
    this.atoms[0] = a;          /* note the anchor atom and the */
    this.bonds[0] = b;          /* (first) bond of the extension */
    this.atoms[1] = d;          /* note the destination atom */
    this.dst = (d.mark >= 0) ? d.mark : this.base.atoms.length;
    if (b.isInRing()            /* if ring extension */
    &&  ((this.mode & RING) != 0)) {
      this.all  = b.getRings(); /* note the ring flags and */
      this.curr = 1;            /* initialize the current flag */
      return this.ring();       /* return whether some ring */
    }                           /* containing this bond is admissible */
    if ((this.mode & BOND) == 0)/* if single bond extension, */
      return false;             /* check for a bond extensions */
    this.atnew = (d.mark < 0) ? 1 : 0;
    this.bdnew = 1;             /* zero/one new atom, one new bond */
    this.size  = 0;             /* clear the extension size */
    this.chns  = this.frag.chns;/* copy the chain counter and */
    return true;                /* return "extension successful" */
  }  /* store() */

  /*------------------------------------------------------------------*/

  protected boolean checkRing ()
  {                             /* --- check a ring extension */
    int  i, n;                  /* loop variable, default index */
    int  sm, dm;                /* indices of source and destination */
    Atom s, a, d;               /* to traverse the ring atoms */
    Bond b, r;                  /* to traverse the ring bonds */

    s = this.atoms[0];          /* get the anchor atom (source), */
    b = this.bonds[0];          /* the first ring bond, */
    d = this.atoms[1];          /* and its destination atom */
    n = this.base.atoms.length; /* get default dest. index */
    for (i = this.size; --i > 0; ) {
      r = this.bonds[i];      /* traverse the other ring bonds that */
      if (r.mark >= 0) return false;         /* are not in the base */
      sm = (r.src.mark >= 0) ? r.src.mark : n;
      dm = (r.dst.mark >= 0) ? r.dst.mark : n;
      if (sm > dm) dm = sm;   /* get index of bond's dest. atom */
      if (dm < this.dst) break;
    }                         /* check for a better first bond */
    if (i > 0) return false;  /* if a better way found, skip ring */
    r = this.bonds[this.size-1];
    if (r.mark >= 0)        return true;
    a  = (r.src  != s) ? r.src  : r.dst;
    dm = (a.mark >= 0) ? a.mark : n;
    if (dm     >  this.dst) return true;
    if (dm     <  this.dst) return false;
    if (r.type >  b.type)   return true;
    if (r.type <  b.type)   return false;
    if (a.type >= d.type)   return true;
    return false;               /* skip rings in wrong direction */
  }  /* checkRing() */

  /*------------------------------------------------------------------*/

  protected void initVars ()
  {                             /* --- init. ring extension variants */
    System.out.println("initVars not implemented in class RightExt");
  }  /* initVars() */

  /*------------------------------------------------------------------*/

  protected boolean variant ()
  {                             /* --- create ring extension variant */
    System.out.println("variant not implemented in class RightExt");
    return true;                /* only print a log message */
  }  /* variant() */

  /*------------------------------------------------------------------*/

  public int compareTo (Fragment frag)
  {                             /* --- compare extension to fragment */
    int       t1, t2;           /* buffers for comparison */
    Embedding emb;              /* to access an embedding */

    if (this.dst < frag.dst) return -1;  /* compare the indices */
    if (this.dst > frag.dst) return +1;  /* of the dest.  atoms */
    if (this.src > frag.src) return -1;  /* compare the indices */
    if (this.src < frag.src) return +1;  /* of the source atoms */
    emb = frag.list;            /* get an underlying embedding */
    t1 = this.bonds[    0   ].type;
    t2 =  emb.bonds[frag.idx].type;
    if (t1 < t2) return -1;     /* compare the types */
    if (t1 > t2) return +1;     /* of the (first) added bond */
    t1 = this.atoms[    1   ].type;
    t2 =  emb.atoms[frag.dst].type;
    if (t1 < t2) return -1;     /* compare the types */
    if (t1 > t2) return +1;     /* of the destination atoms */
    t1 = (this.size < 0) ? 1 : ((this.size > 0) ? -1 : 0);
    t2 = (frag.size < 0) ? 1 : ((frag.size > 0) ? -1 : 0);
    if (t1 < t2) return -1;     /* get the extension types */
    if (t1 > t2) return +1;     /* from the sizes and compare them */
    return (this.size <= 0)     /* compare ring ext. if necessary */
         ? 0 : this.compareRing(frag);
  }  /* compareTo() */

  /*------------------------------------------------------------------*/

  protected boolean unclosable (Fragment frag)
  {                             /* --- check for uncloseable rings */
    int      i, k, n;           /* loop variable, bond counter */
    Molecule mol;               /* the fragment as a molecule */
    Atom     a, d;              /* to traverse the atoms */
    Bond     b;                 /* to traverse the bonds on the path */

    mol = frag.getAsMolecule(); /* get the fragment as a molecule */
    for (i = mol.atomcnt; --i >= 0; )
      mol.atoms[i].mark = i;    /* mark all atoms with their index */
    a = mol.atoms[frag.dst];    /* get rightmost leaf (end of path) */
    while (a.mark > 0) {        /* construct rightmost path upwards */
      for (k = -1, i = a.bondcnt; --i >= 0; ) {
        b = a.bonds[i];         /* traverse the marked bonds */
        d = (b.src != a) ? b.src : b.dst;
        if ((d.mark > k) && (d.mark < a.mark))
          k = d.mark;           /* find the adjacent atom with the */
      }                         /* largest index smaller than own */
      a.mark = -1;              /* mark current atom as extendable */
      a = mol.atoms[k];         /* and get the next path atom */
    }
    a.mark = -1;                /* mark root atom as extendable */
    for (i = mol.atomcnt; --i > 0; ) {
      a = mol.atoms[i];         /* traverse the unextendable atoms */
      if (a.mark < 0) continue; /* (have their indices as markers) */
      a.mark = -1;              /* unmark the atom */
      for (n = 0, k = a.bondcnt; --k >= 0; )
        if ((a.bonds[k].type & Bond.RINGBOND) != 0) n++;
      if (n == 1) return true;  /* if there is a single ring bond, */
    }                           /* a ring cannot be closed anymore */
    return false;               /* all rings may be closable */
  }  /* unclosable() */

  /*------------------------------------------------------------------*/

  protected int isCanonic (int bdi, int ati, int cnt)
  {                             /* --- check prefix words recursively */
    int  i, k, c, m, r;         /* loop variable, atom index, buffer */
    Bond b;                     /* to traverse/access the bonds */
    Atom s, d;                  /* to traverse/access the atoms */

    /* --- check for bonds closing rings --- */