rightext.java

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

    d = this.atoms[cnt-1];      /* get the previous destination atom */
    k = -1; c = (bdi << 2) +1;  /* init. bond and code word index */
    if ((this.word[c]   <  cnt) /* if a new atom is not added (ring) */
    ||  (this.word[c+1] == d.mark)) {  /* or added to rightmost leaf */
      m = (this.word[c] < cnt) ? this.word[c+1] : -1;
      for (i = d.bondcnt; --i >= 0; ) {
        b = d.bonds[i];         /* traverse the unmarked bonds */
        if (b.mark >= 0) continue;    /* of the rightmost leaf */
        s = (b.src != d) ? b.src : b.dst;
        if (s.mark >  m)        /* check for a bond closing a ring */
          return bdi;           /* that should precede the next bond */
        if (s.mark == m) k = i; /* note the index of the bond that */
      }                         /* matches the next reference bond */
    }                           /* (if this bond closes a ring) */

    /* --- process bond closing a ring --- */
    r = this.size;              /* set the default result */
    if (this.word[c] < cnt) {   /* if no new atom is added */
      if (k < 0) return r;      /* check for a corresponding bond */
      b = d.bonds[k];           /* get the corresponding bond */
      if (b.type < this.word[c+2]) return bdi;  /* check the type */
      if (b.type > this.word[c+2]) return r;    /* of the bond */
      if (bdi    >= r-1)           return r;    /* check bond index */
      b.mark =  0;              /* mark the matching bond */
      k = this.isCanonic(bdi+1, 0, cnt);
      b.mark = -1;              /* check remaining bonds recursively, */
      return k;                 /* unmark the matching bond, and */
    }                           /* return the recursion result */

    /* --- process bond adding an atom --- */
    m = this.word[c+1];         /* note index of new source atom */
    while (d.mark > m) {        /* check old path up to new source */
      for (k = -1, i = d.bondcnt; --i >= 0; ) {
        b = d.bonds[i];         /* traverse the marked bonds */
        if (b.mark < 0) return bdi;
        s = (b.src != d) ? b.src : b.dst;
        if ((s.mark > k) && (s.mark < d.mark))
          k = s.mark;           /* find the adjacent atom with the */
      }                         /* largest index smaller than own */
      d = this.atoms[k];        /* get next atom on rightmost path */
    }                           /* (traverse from leaf up to source) */
    s = this.atoms[m];          /* get the current source atom */
    for (i = 0; i < s.bondcnt; i++) {
      b = s.bonds[i];           /* traverse the unmarked bonds */
      if (b.mark >= 0)             continue;
      if (b.type < this.word[c+2]) return bdi;  /* check the type */
      if (b.type > this.word[c+2]) return r;    /* of the bond */
      d = (b.src != s) ? b.src : b.dst;
      if (d.type < this.word[c+3]) return bdi;  /* check the type */
      if (d.type > this.word[c+3]) return r;    /* of the dest. atom */
      if (bdi    >= r-1)           return r;    /* check bond index */
      b.mark = 0;               /* mark the bond and the atom */
      this.atoms[d.mark = cnt] = d;
      k = this.isCanonic(bdi+1, 0, cnt+1);
      b.mark = d.mark = -1;     /* check remaining bonds recursively, */
      if (k < r) { if (k < this.src) return k; r = k; }
    }                           /* evaluate the recursion result */
    return r;                   /* return the overall result */
  }  /* isCanonic() */

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

  protected boolean makeCanonic (int bdi, int ati, int cnt)
  {                             /* --- construct canonic code word */
    int     i, k, c;            /* loop variables, indices */
    Bond    b, x;               /* to traverse/access the bonds */
    Atom    s, d, a;            /* to traverse/access the atoms */
    boolean changed;            /* whether code word was changed */

    c = (bdi << 2) +1;          /* compute index in code word */
    if (bdi >= this.size) {     /* if full code word is constructed */
      i = this.word[c]; this.word[c] = 0;
      return (i != 0);          /* reinstall sentinel and return */
    }                           /* whether code word was changed */

    /* --- bond closing a ring --- */
    b = null; s = null;         /* dummy initialization */
    for ( ; ati < cnt; ati++) { /* search for bonds closing rings */
      d = this.atoms[ati];      /* traverse possible dest. atoms */
      for (i = d.bondcnt; --i >= 0; ) {
        b = d.bonds[i];         /* traverse the unmarked bonds */
        if (b.mark >= 0) continue;
        s = (b.src != d) ? b.src : b.dst;
        if ((s.mark < ati) && (s.mark >= 0))
          break;                /* if there is a backward bond, */
      }                         /* abort the search loop */
      if (i < 0) continue;      /* if no bond found, go to next atom */
      while (--i >= 0) {        /* traverse the remaining unmarked */
        x = d.bonds[i];         /* bonds of the same destination atom */
        if (x.mark >= 0) continue;
        a = (x.src != d) ? x.src : x.dst;
        if ((a.mark < ati) && (a.mark > s.mark)) {
          b = x; s = a; }       /* try to find a source atom */
      }                         /* with a larger index */
      if (d.mark > this.word[c])
        return false;           /* compare destination atom index */
      if (d.mark < this.word[c]) {
        this.word[c]   = d.mark;/* set new destination atom index */
        this.word[c+1] = Integer.MIN_VALUE;
      }                         /* on change invalidate next entry */
      if (s.mark < this.word[c+1])
        return false;           /* compare source atom index */
      if (s.mark > this.word[c+1]) {
        this.word[c+1] = s.mark;/* set new source atom index */
        this.word[c+2] = Integer.MAX_VALUE;
      }                         /* on change invalidate next entry */
      if (b.type > this.word[c+2])
        return false;           /* compare the bond type */
      if (b.type < this.word[c+2]) {
        this.word[c+2] = b.type;/* set new bond type */
        this.word[c+3] = Integer.MAX_VALUE;
      }                         /* on change invalidate next entry */
      if (d.type > this.word[c+3])
        return false;           /* compare destination atom type */
      if (d.type < this.word[c+3]) {
        this.word[c+3] = d.type;/* set new destination atom type */
        this.word[c+4] = Integer.MAX_VALUE;
      }                         /* on change invalidate next entry */
      b.mark =  0;              /* mark the bond and recurse */
      changed = this.makeCanonic(bdi+1, ati, cnt);
      if (changed) this.bonds[bdi] = b;
      b.mark = -1;              /* unmark the bond again */
      return changed;           /* return whether the code word */
    }                           /* was changed in the recursion */

    /* --- bond to a new atom --- */
    if (cnt > this.word[c])     /* if beyond small enough dest., */
      return false;             /* abort the function */
    changed = false;            /* default: code word is unchanged */
    for (k = cnt; --k >= 0; ) { /* traverse possible source atoms */
      s = this.atoms[k];        /* (in reverse order) */
      for (i = s.bondcnt; --i >= 0; )
        if (s.bonds[i].mark < 0) break;
      if (i >= 0) break;        /* if an unmarked bond found, */
    }                           /* abort the search loop */
    if (k < 0) return false;    /* safety check (should not happen) */
    if (cnt < this.word[c]) {   /* if dest. atom has smaller index, */
      this.word[c]   = cnt;     /* replace the code word letter */
      this.word[c+1] = Integer.MIN_VALUE;
    }                           /* invalidate the next entry */
    if (k < this.word[c+1])     /* if source atom has smaller index, */
      return false;             /* the existing code word is better */
    if (k > this.word[c+1]) {   /* if source atom has larger index, */
      this.word[c+1] = k;       /* replace the code word letter */
      this.word[c+2] = Integer.MAX_VALUE;
    }                           /* invalidate the next entry */
    for (i = 0; i < s.bondcnt; i++) {
      b = s.bonds[i];           /* traverse the unmarked bonds */
      if (b.mark >= 0) continue;
      if (b.type > this.word[c+2])
        return changed;         /* compare the bond type */
      if (b.type < this.word[c+2]) {
        this.word[c+2] = b.type;/* set new bond type */
        this.word[c+3] = Integer.MAX_VALUE;
      }                         /* on change invalidate next entry */
      d = (b.src != s) ? b.src : b.dst;
      if (d.type > this.word[c+3])
        return changed;         /* compare destination atom type */
      if (d.type < this.word[c+3]) {
        this.word[c+3] = d.type;/* set new destination atom type */
        this.word[c+4] = Integer.MAX_VALUE;
      }                         /* on change invalidate next entry */
      b.mark = 0;               /* mark the bond and the atom */
      this.atoms[d.mark = cnt] = d;
      if (this.makeCanonic(bdi+1, ati, cnt+1)) {
        this.bonds[bdi] = b; changed = true; }
      b.mark = d.mark = -1;     /* recursively construct code word, */
    }                           /* then unmark bond and atom again */
    return changed;             /* return whether bonds were replaced */
  }  /* makeCanonic() */

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

  protected boolean adaptRing (Fragment frag)
  {                             /* --- adapt and check ring extension */
    System.out.println("adaptRings not implemented in class RightExt");
    return true;                /* only print a log message */
  }  /* adaptRing() */

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

  protected int makeWord (Molecule mol, int bondcnt)
  {                             /* --- construct (part of) code word */
    int   i, k, n;              /* loop variables, buffers */
    Bond  b;                    /* to traverse the bonds */

    n = (bondcnt << 2) +1;      /* compute number of characters */ 
    if (n > this.word.length)   /* if the word gets too long, */
      this.word = new int[n];   /* create a new buffer */
    for (i = mol.atomcnt; --i >= 0; )
      mol.atoms[i].mark = i;    /* number the atoms of the molecule */
    this.word[k = 0] = mol.atoms[0].type;
    for (i = 0; i < bondcnt; i++) {
      b = mol.bonds[i];         /* traverse the molecule's bonds */
      if (b.src.mark < b.dst.mark) {
        this.word[++k] = b.dst.mark; /* if "forward" bond */
        this.word[++k] = b.src.mark; /* start with source atom */
        this.word[++k] = b.type;
        this.word[++k] = b.dst.type; }
      else {                         /* if "backward" bond */
        this.word[++k] = b.src.mark; /* start with destination atom */
        this.word[++k] = b.dst.mark;
        this.word[++k] = b.type;
        this.word[++k] = b.src.type;
      }                         /* describe a bond of the molecule */
    }
    return n;                   /* return the number of characters */
  }  /* makeWord() */

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

  protected int compareWord (Molecule mol, int bondcnt)
  {                             /* --- compare molecule to code word */
    int  i, k;                  /* loop variables */
    Bond b;                     /* to traverse the bonds */
    Atom s, d;                  /* to traverse the atoms */

    for (i = mol.atomcnt; --i >= 0; )
      mol.atoms[i].mark = i;    /* number the atoms of the molecule */
    s = mol.atoms[0];           /* compare the type of the root atom */
    if (s.type < this.word[0]) return -1;
    if (s.type > this.word[0]) return +1;
    for (i = k = 0; i < bondcnt; i++) {
      b = mol.bonds[i];         /* traverse the molecule's bonds */
      s = b.src; d = b.dst;     /* get the source and dest. atoms */
      if (s.mark > d.mark) { s = d; d = b.src; }
      if (d.mark < this.word[++k]) return -1;
      if (d.mark > this.word[  k]) return +1;
      if (s.mark > this.word[++k]) return -1;
      if (s.mark < this.word[  k]) return +1;
      if (b.type < this.word[++k]) return -1;
      if (b.type > this.word[  k]) return +1;
      if (d.type < this.word[++k]) return -1;
      if (d.type > this.word[  k]) return +1;
    }                           /* return sign of difference */
    return 0;                   /* otherwise return 'equal' */
  }  /* compareWord() */

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

  protected String makeString (Molecule mol, boolean create)
  {                             /* --- print a molecule's code word */
    int          i, k, n;       /* loop variable, buffers */
    StringBuffer s;             /* created description */

    if (create)                 /* construct molecule's code word */
      this.makeWord(mol, mol.bondcnt);
    k = this.word[0];           /* get and decode type of root atom */
    if (mol.map != null) k = mol.map.decode(k);
    s = new StringBuffer(Atom.names[k & Atom.TYPEMASK]);
    n = mol.bondcnt << 2;       /* get the number of characters */
    for (i = 0; i < n; ) {      /* traverse the characters */
      s.append("|" +this.word[++i]);   /* destination atom index */
      s.append(" " +this.word[++i]);   /* source      atom index */
      s.append(Bond.names[this.word[++i] & Bond.TYPEMASK]);
      k = this.word[++i];       /* bond and atom type */
      if (mol.map != null) k = mol.map.decode(k);
      s.append(Atom.names[k & Atom.TYPEMASK]);
    }                           /* store the edge descriptions */
    return s.toString();        /* return created string description */
  }  /* makeString() */

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

  protected int compareBond (Bond b1, Bond b2, int next)
  {                             /* --- compare two bonds */
    Atom s1, d1, s2, d2;        /* source and destination atoms */
    int  i1, i2;                /* indices of the destination atoms */

    s1 = b1.src; d1 = b1.dst;   /* get atoms of first bond */
    if ((d1.mark >= 0) && (d1.mark < s1.mark)) {
      s1 = d1; d1 = b1.src; }   /* exchange atoms if necessary */
    s2 = b2.src; d2 = b2.dst;   /* get atoms of second bond */
    if ((d2.mark >= 0) && (d2.mark < s2.mark)) {
      s2 = d2; d2 = b2.src; }   /* exchange atoms if necessary */
    i1 = (d1.mark >= 0) ? d1.mark : next;
    i2 = (d2.mark >= 0) ? d2.mark : next;
    if (i1      < i2)      return -1;  /* compare the indices */
    if (i1      > i2)      return +1;  /* of the destination atoms */
    if (s1.mark > s2.mark) return -1;  /* compare the indices */
    if (s1.mark < s2.mark) return +1;  /* of the source atoms */
    if (b1.type < b2.type) return -1;  /* compare the types */
    if (b1.type > b2.type) return +1;  /* of the two bonds */
    if (d1.type < d2.type) return -1;  /* compare the types */ 
    if (d1.type > d2.type) return +1;  /* of the destination atoms */
    return 0;                   /* otherwise the bonds are equal */
  }  /* compareBond() */

}  /* class RightExt */