removing cyclic call
[jalview.git] / src / jalview / datamodel / Sequence.java
index 783b1bc..6e9e1cd 100755 (executable)
 package jalview.datamodel;
 
 import jalview.analysis.AlignSeq;
-import jalview.api.DBRefEntryI;
 import jalview.datamodel.features.SequenceFeatures;
 import jalview.datamodel.features.SequenceFeaturesI;
 import jalview.util.Comparison;
 import jalview.util.DBRefUtils;
 import jalview.util.MapList;
 import jalview.util.StringUtils;
+import jalview.ws.params.InvalidArgumentException;
 
 import java.util.ArrayList;
 import java.util.Arrays;
 import java.util.BitSet;
 import java.util.Collections;
 import java.util.Enumeration;
+import java.util.Iterator;
 import java.util.List;
 import java.util.ListIterator;
 import java.util.Vector;
 
-import com.stevesoft.pat.Regex;
-
 import fr.orsay.lri.varna.models.rna.RNA;
 
 /**
  * 
- * Implements the SequenceI interface for a char[] based sequence object.
- * 
- * @author $author$
- * @version $Revision$
+ * Implements the SequenceI interface for a char[] based sequence object
  */
 public class Sequence extends ASequence implements SequenceI
 {
-  private static final Regex limitrx = new Regex(
-          "[/][0-9]{1,}[-][0-9]{1,}$");
 
-  private static final Regex endrx = new Regex("[0-9]{1,}$");
+       /**
+        * A subclass that gives us access to modCount, which tracks 
+        * whether there have been any changes. We use this to update 
+        * @author hansonr
+        *
+        * @param <T>
+        */
+  @SuppressWarnings("serial")
+  protected class DBModList<T> extends ArrayList<DBRefEntry> {
+
+         protected int getModCount() {
+                 return modCount;
+         }
+         
+  }
 
-  SequenceI datasetSequence;
+SequenceI datasetSequence;
 
   String name;
 
@@ -72,10 +80,17 @@ public class Sequence extends ASequence implements SequenceI
 
   String vamsasId;
 
-  DBRefEntry[] dbrefs;
+  private DBModList<DBRefEntry> dbrefs; // controlled acces
 
-  RNA rna;
+  /**
+   * a flag to let us know that elements have changed in dbrefs
+   * 
+   * @author Bob Hanson
+   */
+  private int refModCount = 0;
 
+  RNA rna;
+  
   /**
    * This annotation is displayed below the alignment but the positions are tied
    * to the residues of this sequence
@@ -84,12 +99,7 @@ public class Sequence extends ASequence implements SequenceI
    */
   Vector<AlignmentAnnotation> annotation;
 
-  /**
-   * The index of the sequence in a MSA
-   */
-  int index = -1;
-
-  private SequenceFeatures sequenceFeatureStore;
+  private SequenceFeaturesI sequenceFeatureStore;
 
   /*
    * A cursor holding the approximate current view position to the sequence,
@@ -151,25 +161,49 @@ public class Sequence extends ASequence implements SequenceI
     checkValidRange();
   }
 
+  /**
+   * If 'name' ends in /i-j, where i >= j > 0 are integers, extracts i and j as
+   * start and end respectively and removes the suffix from the name
+   */
   void parseId()
   {
     if (name == null)
     {
-      System.err
-              .println("POSSIBLE IMPLEMENTATION ERROR: null sequence name passed to constructor.");
+      System.err.println(
+              "POSSIBLE IMPLEMENTATION ERROR: null sequence name passed to constructor.");
       name = "";
     }
-    // Does sequence have the /start-end signature?
-    if (limitrx.search(name))
+    int slashPos = name.lastIndexOf('/');
+    if (slashPos > -1 && slashPos < name.length() - 1)
     {
-      name = limitrx.left();
-      endrx.search(limitrx.stringMatched());
-      setStart(Integer.parseInt(limitrx.stringMatched().substring(1,
-              endrx.matchedFrom() - 1)));
-      setEnd(Integer.parseInt(endrx.stringMatched()));
+      String suffix = name.substring(slashPos + 1);
+      String[] range = suffix.split("-");
+      if (range.length == 2)
+      {
+        try
+        {
+          int from = Integer.valueOf(range[0]);
+          int to = Integer.valueOf(range[1]);
+          if (from > 0 && to >= from)
+          {
+            name = name.substring(0, slashPos);
+            setStart(from);
+            setEnd(to);
+            checkValidRange();
+          }
+        } catch (NumberFormatException e)
+        {
+          // leave name unchanged if suffix is invalid
+        }
+      }
     }
   }
 
+  /**
+   * Ensures that 'end' is not before the end of the sequence, that is,
+   * (end-start+1) is at least as long as the count of ungapped positions. Note
+   * that end is permitted to be beyond the end of the sequence data.
+   */
   void checkValidRange()
   {
     // Note: JAL-774 :
@@ -178,7 +212,7 @@ public class Sequence extends ASequence implements SequenceI
       int endRes = 0;
       for (int j = 0; j < sequence.length; j++)
       {
-        if (!jalview.util.Comparison.isGap(sequence[j]))
+        if (!Comparison.isGap(sequence[j]))
         {
           endRes++;
         }
@@ -274,12 +308,12 @@ public class Sequence extends ASequence implements SequenceI
      */
     if (datasetSequence == null)
     {
-      if (seq.getDBRefs() != null)
+      List<DBRefEntry> dbr = seq.getDBRefs();
+      if (dbr != null)
       {
-        DBRefEntry[] dbr = seq.getDBRefs();
-        for (int i = 0; i < dbr.length; i++)
+        for (int i = 0, n = dbr.size(); i < n; i++)
         {
-          addDBRef(new DBRefEntry(dbr[i]));
+          addDBRef(new DBRefEntry(dbr.get(i)));
         }
       }
 
@@ -396,7 +430,7 @@ public class Sequence extends ASequence implements SequenceI
   {
     if (pdbIds == null)
     {
-      pdbIds = new Vector<PDBEntry>();
+      pdbIds = new Vector<>();
       pdbIds.add(entry);
       return true;
     }
@@ -432,7 +466,7 @@ public class Sequence extends ASequence implements SequenceI
   @Override
   public Vector<PDBEntry> getAllPDBEntries()
   {
-    return pdbIds == null ? new Vector<PDBEntry>() : pdbIds;
+    return pdbIds == null ? new Vector<>() : pdbIds;
   }
 
   /**
@@ -453,15 +487,15 @@ public class Sequence extends ASequence implements SequenceI
   }
 
   /**
-   * DOCUMENT ME!
+   * Sets the sequence name. If the name ends in /start-end, then the start-end
+   * values are parsed out and set, and the suffix is removed from the name.
    * 
-   * @param name
-   *          DOCUMENT ME!
+   * @param theName
    */
   @Override
-  public void setName(String name)
+  public void setName(String theName)
   {
-    this.name = name;
+    this.name = theName;
     this.parseId();
   }
 
@@ -645,10 +679,10 @@ public class Sequence extends ASequence implements SequenceI
   }
 
   /**
-   * DOCUMENT ME!
+   * Sets the sequence description, and also parses out any special formats of
+   * interest
    * 
    * @param desc
-   *          DOCUMENT ME!
    */
   @Override
   public void setDescription(String desc)
@@ -656,10 +690,67 @@ public class Sequence extends ASequence implements SequenceI
     this.description = desc;
   }
 
+  @Override
+  public void setGeneLoci(String speciesId, String assemblyId,
+          String chromosomeId, MapList map)
+  {
+    addDBRef(new DBRefEntry(speciesId, assemblyId, DBRefEntry.CHROMOSOME
+            + ":" + chromosomeId, new Mapping(map)));
+  }
+
   /**
-   * DOCUMENT ME!
+   * Returns the gene loci mapping for the sequence (may be null)
    * 
-   * @return DOCUMENT ME!
+   * @return
+   */
+  @Override
+  public GeneLociI getGeneLoci()
+  {
+    List<DBRefEntry> refs = getDBRefs();
+    if (refs != null)
+    {
+      for (final DBRefEntry ref : refs)
+      {
+        if (ref.isChromosome())
+        {
+          return new GeneLociI()
+          {
+            @Override
+            public String getSpeciesId()
+            {
+              return ref.getSource();
+            }
+
+            @Override
+            public String getAssemblyId()
+            {
+              return ref.getVersion();
+            }
+
+            @Override
+            public String getChromosomeId()
+            {
+              // strip off "chromosome:" prefix to chrId
+              return ref.getAccessionId().substring(
+                      DBRefEntry.CHROMOSOME.length() + 1);
+            }
+
+            @Override
+            public MapList getMap()
+            {
+              return ref.getMap().getMap();
+            }
+          };
+        }
+      }
+    }
+    return null;
+  }
+
+  /**
+   * Answers the description
+   * 
+   * @return
    */
   @Override
   public String getDescription()
@@ -727,6 +818,7 @@ public class Sequence extends ASequence implements SequenceI
      * preserve end residue column provided cursor was valid
      */
     int endColumn = isValidCursor(cursor) ? cursor.lastColumnPosition : 0;
+
     if (residuePos == this.end)
     {
       endColumn = column;
@@ -745,7 +837,7 @@ public class Sequence extends ASequence implements SequenceI
    * @param curs
    * @return
    */
-  protected int findIndex(int pos, SequenceCursor curs)
+  protected int findIndex(final int pos, SequenceCursor curs)
   {
     if (!isValidCursor(curs))
     {
@@ -763,16 +855,20 @@ public class Sequence extends ASequence implements SequenceI
     /*
      * move left or right to find pos from hint.position
      */
-    int col = curs.columnPosition - 1; // convert from base 1 to 0-based array
-                                       // index
+    int col = curs.columnPosition - 1; // convert from base 1 to base 0
     int newPos = curs.residuePosition;
     int delta = newPos > pos ? -1 : 1;
 
     while (newPos != pos)
     {
       col += delta; // shift one column left or right
-      if (col < 0 || col == sequence.length)
+      if (col < 0)
+      {
+        break;
+      }
+      if (col == sequence.length)
       {
+        col--; // return last column if we failed to reach pos
         break;
       }
       if (!Comparison.isGap(sequence[col]))
@@ -782,7 +878,14 @@ public class Sequence extends ASequence implements SequenceI
     }
 
     col++; // convert back to base 1
-    updateCursor(pos, col, curs.firstColumnPosition);
+
+    /*
+     * only update cursor if we found the target position
+     */
+    if (newPos == pos)
+    {
+      updateCursor(pos, col, curs.firstColumnPosition);
+    }
 
     return col;
   }
@@ -984,6 +1087,53 @@ public class Sequence extends ASequence implements SequenceI
   }
 
   /**
+   * {@inheritDoc}
+   */
+  @Override
+  public Range findPositions(int fromColumn, int toColumn)
+  {
+    if (toColumn < fromColumn || fromColumn < 1)
+    {
+      return null;
+    }
+
+    /*
+     * find the first non-gapped position, if any
+     */
+    int firstPosition = 0;
+    int col = fromColumn - 1;
+    int length = sequence.length;
+    while (col < length && col < toColumn)
+    {
+      if (!Comparison.isGap(sequence[col]))
+      {
+        firstPosition = findPosition(col++);
+        break;
+      }
+      col++;
+    }
+
+    if (firstPosition == 0)
+    {
+      return null;
+    }
+
+    /*
+     * find the last non-gapped position
+     */
+    int lastPosition = firstPosition;
+    while (col < length && col < toColumn)
+    {
+      if (!Comparison.isGap(sequence[col++]))
+      {
+        lastPosition++;
+      }
+    }
+
+    return new Range(firstPosition, lastPosition);
+  }
+
+  /**
    * Returns an int array where indices correspond to each residue in the
    * sequence and the element value gives its position in the alignment
    * 
@@ -1012,6 +1162,27 @@ public class Sequence extends ASequence implements SequenceI
     return map;
   }
 
+  /**
+   * Build a bitset corresponding to sequence gaps
+   * 
+   * @return a BitSet where set values correspond to gaps in the sequence
+   */
+  @Override
+  public BitSet gapBitset()
+  {
+    BitSet gaps = new BitSet(sequence.length);
+    int j = 0;
+    while (j < sequence.length)
+    {
+      if (jalview.util.Comparison.isGap(sequence[j]))
+      {
+        gaps.set(j);
+      }
+      j++;
+    }
+    return gaps;
+  }
+
   @Override
   public int[] findPositionMap()
   {
@@ -1035,7 +1206,7 @@ public class Sequence extends ASequence implements SequenceI
   @Override
   public List<int[]> getInsertions()
   {
-    ArrayList<int[]> map = new ArrayList<int[]>();
+    ArrayList<int[]> map = new ArrayList<>();
     int lastj = -1, j = 0;
     int pos = start;
     int seqlen = sequence.length;
@@ -1101,7 +1272,7 @@ public class Sequence extends ASequence implements SequenceI
   }
 
   @Override
-  public void deleteChars(int i, int j)
+  public void deleteChars(final int i, final int j)
   {
     int newstart = start, newend = end;
     if (i >= sequence.length || i < 0)
@@ -1113,62 +1284,76 @@ public class Sequence extends ASequence implements SequenceI
     boolean createNewDs = false;
     // TODO: take a (second look) at the dataset creation validation method for
     // the very large sequence case
-    int eindex = -1, sindex = -1;
-    boolean ecalc = false, scalc = false;
-    for (int s = i; s < j; s++)
+
+    int startIndex = findIndex(start) - 1;
+    int endIndex = findIndex(end) - 1;
+    int startDeleteColumn = -1; // for dataset sequence deletions
+    int deleteCount = 0;
+
+    for (int s = i; s < j && s < sequence.length; s++)
     {
-      if (jalview.schemes.ResidueProperties.aaIndex[sequence[s]] != 23)
+      if (Comparison.isGap(sequence[s]))
+      {
+        continue;
+      }
+      deleteCount++;
+      if (startDeleteColumn == -1)
+      {
+        startDeleteColumn = findPosition(s) - start;
+      }
+      if (createNewDs)
       {
-        if (createNewDs)
+        newend--;
+      }
+      else
+      {
+        if (startIndex == s)
         {
-          newend--;
+          /*
+           * deleting characters from start of sequence; new start is the
+           * sequence position of the next column (position to the right
+           * if the column position is gapped)
+           */
+          newstart = findPosition(j);
+          break;
         }
         else
         {
-          if (!scalc)
+          if (endIndex < j)
           {
-            sindex = findIndex(start) - 1;
-            scalc = true;
-          }
-          if (sindex == s)
-          {
-            // delete characters including start of sequence
-            newstart = findPosition(j);
-            break; // don't need to search for any more residue characters.
+            /*
+             * deleting characters at end of sequence; new end is the sequence
+             * position of the column before the deletion; subtract 1 if this is
+             * gapped since findPosition returns the next sequence position
+             */
+            newend = findPosition(i - 1);
+            if (Comparison.isGap(sequence[i - 1]))
+            {
+              newend--;
+            }
+            break;
           }
           else
           {
-            // delete characters after start.
-            if (!ecalc)
-            {
-              eindex = findIndex(end) - 1;
-              ecalc = true;
-            }
-            if (eindex < j)
-            {
-              // delete characters at end of sequence
-              newend = findPosition(i - 1);
-              break; // don't need to search for any more residue characters.
-            }
-            else
-            {
-              createNewDs = true;
-              newend--; // decrease end position by one for the deleted residue
-              // and search further
-            }
+            createNewDs = true;
+            newend--;
           }
         }
       }
     }
-    // deletion occured in the middle of the sequence
+
     if (createNewDs && this.datasetSequence != null)
     {
-      // construct a new sequence
+      /*
+       * if deletion occured in the middle of the sequence,
+       * construct a new dataset sequence and delete the residues
+       * that were deleted from the aligned sequence
+       */
       Sequence ds = new Sequence(datasetSequence);
+      ds.deleteChars(startDeleteColumn, startDeleteColumn + deleteCount);
+      datasetSequence = ds;
       // TODO: remove any non-inheritable properties ?
       // TODO: create a sequence mapping (since there is a relation here ?)
-      ds.deleteChars(i, j);
-      datasetSequence = ds;
     }
     start = newstart;
     end = newend;
@@ -1225,24 +1410,35 @@ public class Sequence extends ASequence implements SequenceI
     vamsasId = id;
   }
 
-  @Override
-  public void setDBRefs(DBRefEntry[] dbref)
+  @SuppressWarnings("deprecation")
+@Override
+  public void setDBRefs(List<DBRefEntry> newDBrefs) throws InvalidArgumentException
   {
     if (dbrefs == null && datasetSequence != null
             && this != datasetSequence)
     {
-      datasetSequence.setDBRefs(dbref);
+      datasetSequence.setDBRefs((DBModList<DBRefEntry>)newDBrefs);
       return;
     }
-    dbrefs = dbref;
-    if (dbrefs != null)
-    {
-      DBRefUtils.ensurePrimaries(this);
-    }
+    if (newDBrefs != null && !(newDBrefs instanceof DBModList<?>))
+       throw new InvalidArgumentException("DBrefs must have DBModList class");
+
+    dbrefs = (DBModList<DBRefEntry>)newDBrefs;
+    refModCount = 0;
+  }
+
+  @Override
+  public void getDBRefsFrom(SequenceI seq) {
+         try {
+               setDBRefs(seq.getDBRefs());
+       } catch (InvalidArgumentException e) {
+               // TODO Auto-generated catch block
+               e.printStackTrace();
+       }
   }
 
   @Override
-  public DBRefEntry[] getDBRefs()
+  public List<DBRefEntry> getDBRefs()
   {
     if (dbrefs == null && datasetSequence != null
             && this != datasetSequence)
@@ -1252,6 +1448,7 @@ public class Sequence extends ASequence implements SequenceI
     return dbrefs;
   }
 
+
   @Override
   public void addDBRef(DBRefEntry entry)
   {
@@ -1263,12 +1460,12 @@ public class Sequence extends ASequence implements SequenceI
 
     if (dbrefs == null)
     {
-      dbrefs = new DBRefEntry[0];
+      dbrefs = new DBModList<DBRefEntry>();    
     }
 
-    for (DBRefEntryI dbr : dbrefs)
+    for (int ib = 0, nb= dbrefs.size(); ib < nb; ib++)
     {
-      if (dbr.updateFrom(entry))
+      if (dbrefs.get(ib).updateFrom(entry))
       {
         /*
          * found a dbref that either matched, or could be
@@ -1278,18 +1475,20 @@ public class Sequence extends ASequence implements SequenceI
       }
     }
 
-    /*
-     * extend the array to make room for one more
-     */
-    // TODO use an ArrayList instead
-    int j = dbrefs.length;
-    DBRefEntry[] temp = new DBRefEntry[j + 1];
-    System.arraycopy(dbrefs, 0, temp, 0, j);
-    temp[temp.length - 1] = entry;
 
-    dbrefs = temp;
-
-    DBRefUtils.ensurePrimaries(this);
+//    ///  BH OUCH!
+//    /*
+//     * extend the array to make room for one more
+//     */
+//    // TODO use an ArrayList instead
+//    int j = dbrefs.length;
+//    List<DBRefEntry> temp = new DBRefEntry[j + 1];
+//    System.arraycopy(dbrefs, 0, temp, 0, j);
+//    temp[temp.length - 1] = entry;
+//
+//    dbrefs = temp;
+    
+    dbrefs.add(entry);
   }
 
   @Override
@@ -1317,8 +1516,9 @@ public class Sequence extends ASequence implements SequenceI
   @Override
   public AlignmentAnnotation[] getAnnotation()
   {
-    return annotation == null ? null : annotation
-            .toArray(new AlignmentAnnotation[annotation.size()]);
+    return annotation == null ? null
+            : annotation
+                    .toArray(new AlignmentAnnotation[annotation.size()]);
   }
 
   @Override
@@ -1332,7 +1532,7 @@ public class Sequence extends ASequence implements SequenceI
   {
     if (this.annotation == null)
     {
-      this.annotation = new Vector<AlignmentAnnotation>();
+      this.annotation = new Vector<>();
     }
     if (!this.annotation.contains(annotation))
     {
@@ -1401,6 +1601,8 @@ public class Sequence extends ASequence implements SequenceI
 
   private int _seqhash = 0;
 
+private List<DBRefEntry> primaryRefs;
+
   /**
    * Answers false if the sequence is more than 85% nucleotide (ACGTU), else
    * true
@@ -1430,8 +1632,9 @@ public class Sequence extends ASequence implements SequenceI
   {
     if (datasetSequence == null)
     {
-      Sequence dsseq = new Sequence(getName(), AlignSeq.extractGaps(
-              jalview.util.Comparison.GapChars, getSequenceAsString()),
+      Sequence dsseq = new Sequence(getName(),
+              AlignSeq.extractGaps(jalview.util.Comparison.GapChars,
+                      getSequenceAsString()),
               getStart(), getEnd());
 
       datasetSequence = dsseq;
@@ -1498,7 +1701,7 @@ public class Sequence extends ASequence implements SequenceI
       return null;
     }
 
-    Vector<AlignmentAnnotation> subset = new Vector<AlignmentAnnotation>();
+    Vector<AlignmentAnnotation> subset = new Vector<>();
     Enumeration<AlignmentAnnotation> e = annotation.elements();
     while (e.hasMoreElements())
     {
@@ -1531,13 +1734,14 @@ public class Sequence extends ASequence implements SequenceI
       // TODO: could merge DBRefs
       return datasetSequence.updatePDBIds();
     }
-    if (dbrefs == null || dbrefs.length == 0)
+    if (dbrefs == null || dbrefs.size() == 0)
     {
       return false;
     }
     boolean added = false;
-    for (DBRefEntry dbr : dbrefs)
+    for (int ib = 0, nb = dbrefs.size(); ib < nb; ib++)
     {
+      DBRefEntry dbr = dbrefs.get(ib);
       if (DBRefSource.PDB.equals(dbr.getSource()))
       {
         /*
@@ -1573,7 +1777,7 @@ public class Sequence extends ASequence implements SequenceI
       List<SequenceFeature> sfs = entry.getSequenceFeatures();
       for (SequenceFeature feature : sfs)
       {
-        SequenceFeature sf[] = (mp != null) ? mp.locateFeature(feature)
+       SequenceFeature sf[] = (mp != null) ? mp.locateFeature(feature)
                 : new SequenceFeature[] { new SequenceFeature(feature) };
         if (sf != null)
         {
@@ -1596,12 +1800,12 @@ public class Sequence extends ASequence implements SequenceI
       }
     }
     // transfer database references
-    DBRefEntry[] entryRefs = entry.getDBRefs();
+    List<DBRefEntry> entryRefs = entry.getDBRefs();
     if (entryRefs != null)
     {
-      for (int r = 0; r < entryRefs.length; r++)
+      for (int r = 0, n = entryRefs.size(); r < n; r++)
       {
-        DBRefEntry newref = new DBRefEntry(entryRefs[r]);
+        DBRefEntry newref = new DBRefEntry(entryRefs.get(r));
         if (newref.getMap() != null && mp != null)
         {
           // remap ref using our local mapping
@@ -1616,30 +1820,6 @@ public class Sequence extends ASequence implements SequenceI
     }
   }
 
-  /**
-   * @return The index (zero-based) on this sequence in the MSA. It returns
-   *         {@code -1} if this information is not available.
-   */
-  @Override
-  public int getIndex()
-  {
-    return index;
-  }
-
-  /**
-   * Defines the position of this sequence in the MSA. Use the value {@code -1}
-   * if this information is undefined.
-   * 
-   * @param The
-   *          position for this sequence. This value is zero-based (zero for
-   *          this first sequence)
-   */
-  @Override
-  public void setIndex(int value)
-  {
-    index = value;
-  }
-
   @Override
   public void setRNA(RNA r)
   {
@@ -1656,7 +1836,7 @@ public class Sequence extends ASequence implements SequenceI
   public List<AlignmentAnnotation> getAlignmentAnnotations(String calcId,
           String label)
   {
-    List<AlignmentAnnotation> result = new ArrayList<AlignmentAnnotation>();
+    List<AlignmentAnnotation> result = new ArrayList<>();
     if (this.annotation != null)
     {
       for (AlignmentAnnotation ann : annotation)
@@ -1699,6 +1879,8 @@ public class Sequence extends ASequence implements SequenceI
     return null;
   }
 
+  private List<DBRefEntry> tmpList;
+  
   @Override
   public List<DBRefEntry> getPrimaryDBRefs()
   {
@@ -1706,16 +1888,24 @@ public class Sequence extends ASequence implements SequenceI
     {
       return datasetSequence.getPrimaryDBRefs();
     }
-    if (dbrefs == null || dbrefs.length == 0)
+    if (dbrefs == null || dbrefs.size() == 0)
     {
       return Collections.emptyList();
     }
     synchronized (dbrefs)
     {
-      List<DBRefEntry> primaries = new ArrayList<DBRefEntry>();
-      DBRefEntry[] tmp = new DBRefEntry[1];
-      for (DBRefEntry ref : dbrefs)
+      if (refModCount == dbrefs.getModCount() && primaryRefs != null)
+         return primaryRefs; // no changes
+      refModCount = dbrefs.getModCount(); 
+      List<DBRefEntry> primaries = (primaryRefs == null ? (primaryRefs = new ArrayList<>()) : primaryRefs);
+      primaries.clear();
+      if (tmpList == null) {
+         tmpList = new ArrayList<>();
+         tmpList.add(null); // for replacement 
+      }
+      for (int i = 0, n = dbrefs.size(); i < n; i++)
       {
+       DBRefEntry ref = dbrefs.get(i);
         if (!ref.isPrimaryCandidate())
         {
           continue;
@@ -1730,8 +1920,7 @@ public class Sequence extends ASequence implements SequenceI
           }
         }
         // whilst it looks like it is a primary ref, we also sanity check type
-        if (DBRefUtils.getCanonicalName(DBRefSource.PDB).equals(
-                DBRefUtils.getCanonicalName(ref.getSource())))
+        if (DBRefSource.PDB_CANONICAL_NAME.equals(ref.getCanonicalSourceName()))
         {
           // PDB dbrefs imply there should be a PDBEntry associated
           // TODO: tighten PDB dbrefs
@@ -1740,21 +1929,23 @@ public class Sequence extends ASequence implements SequenceI
           // handle on the PDBEntry, and a real mapping between sequence and
           // extracted sequence from PDB file
           PDBEntry pdbentry = getPDBEntry(ref.getAccessionId());
-          if (pdbentry != null && pdbentry.getFile() != null)
+          if (pdbentry == null || pdbentry.getFile() == null)
           {
-            primaries.add(ref);
+            continue;
           }
-          continue;
-        }
-        // check standard protein or dna sources
-        tmp[0] = ref;
-        DBRefEntry[] res = DBRefUtils.selectDbRefs(!isProtein(), tmp);
-        if (res != null && res[0] == tmp[0])
-        {
-          primaries.add(ref);
-          continue;
-        }
+        } else {
+             // check standard protein or dna sources
+             tmpList.set(0, ref);
+             List<DBRefEntry> res = DBRefUtils.selectDbRefs(!isProtein(), tmpList);
+             if (res == null || res.get(0) != tmpList.get(0))
+             {
+               continue;
+             }
+           }
+        primaries.add(ref);
       }
+      
+      DBRefUtils.ensurePrimaries(this, primaries);
       return primaries;
     }
   }
@@ -1767,14 +1958,11 @@ public class Sequence extends ASequence implements SequenceI
           String... types)
   {
     int startPos = findPosition(fromColumn - 1); // convert base 1 to base 0
-    int endPos = findPosition(toColumn - 1);
-    // to trace / debug behaviour:
-    // System.out
-    // .println(String
-    // .format("%s.findFeatures columns [%d-%d] positions [%d-%d] leaves cursor %s",
-    // getName(), fromColumn, toColumn, startPos,
-    // endPos, cursor));
-    List<SequenceFeature> result = new ArrayList<>();
+    int endPos = fromColumn == toColumn ? startPos
+            : findPosition(toColumn - 1);
+
+    List<SequenceFeature> result = getFeatures().findFeatures(startPos,
+            endPos, types);
     if (datasetSequence != null)
     {
       result = datasetSequence.getFeatures().findFeatures(startPos, endPos,
@@ -1786,33 +1974,40 @@ public class Sequence extends ASequence implements SequenceI
     }
 
     /*
-     * if the start or end column is gapped, startPos or endPos may be to the 
-     * left or right, and we may have included adjacent or enclosing features;
+     * if end column is gapped, endPos may be to the right, 
+     * and we may have included adjacent or enclosing features;
      * remove any that are not enclosing, non-contact features
      */
-    if (endPos > this.end || Comparison.isGap(sequence[fromColumn - 1])
-            || Comparison.isGap(sequence[toColumn - 1]))
+    boolean endColumnIsGapped = toColumn > 0 && toColumn <= sequence.length
+            && Comparison.isGap(sequence[toColumn - 1]);
+    if (endPos > this.end || endColumnIsGapped)
     {
       ListIterator<SequenceFeature> it = result.listIterator();
       while (it.hasNext())
       {
         SequenceFeature sf = it.next();
-        int featureStartColumn = findIndex(sf.getBegin());
-        int featureEndColumn = findIndex(sf.getEnd());
-        boolean noOverlap = featureStartColumn > toColumn
-                        || featureEndColumn < fromColumn;
-
-        /*
-         * reject an 'enclosing' feature if it is actually a contact feature
-         */
-        if (sf.isContactFeature() && featureStartColumn < fromColumn
-                && featureEndColumn > toColumn)
+        int sfBegin = sf.getBegin();
+        int sfEnd = sf.getEnd();
+        int featureStartColumn = findIndex(sfBegin);
+        if (featureStartColumn > toColumn)
         {
-          noOverlap = true;
+          it.remove();
         }
-        if (noOverlap)
+        else if (featureStartColumn < fromColumn)
         {
-          it.remove();
+          int featureEndColumn = sfEnd == sfBegin ? featureStartColumn
+                  : findIndex(sfEnd);
+          if (featureEndColumn < fromColumn)
+          {
+            it.remove();
+          }
+          else if (featureEndColumn > toColumn && sf.isContactFeature())
+          {
+            /*
+             * remove an enclosing feature if it is a contact feature
+             */
+            it.remove();
+          }
         }
       }
     }
@@ -1859,4 +2054,73 @@ public class Sequence extends ASequence implements SequenceI
 
     return count;
   }
+
+  @Override
+  public String getSequenceStringFromIterator(Iterator<int[]> it)
+  {
+    StringBuilder newSequence = new StringBuilder();
+    while (it.hasNext())
+    {
+      int[] block = it.next();
+      if (it.hasNext())
+      {
+        newSequence.append(getSequence(block[0], block[1] + 1));
+      }
+      else
+      {
+        newSequence.append(getSequence(block[0], block[1]));
+      }
+    }
+
+    return newSequence.toString();
+  }
+
+  @Override
+  public int firstResidueOutsideIterator(Iterator<int[]> regions)
+  {
+    int start = 0;
+
+    if (!regions.hasNext())
+    {
+      return findIndex(getStart()) - 1;
+    }
+
+    // Simply walk along the sequence whilst watching for region
+    // boundaries
+    int hideStart = getLength();
+    int hideEnd = -1;
+    boolean foundStart = false;
+
+    // step through the non-gapped positions of the sequence
+    for (int i = getStart(); i <= getEnd() && (!foundStart); i++)
+    {
+      // get alignment position of this residue in the sequence
+      int p = findIndex(i) - 1;
+
+      // update region start/end
+      while (hideEnd < p && regions.hasNext())
+      {
+        int[] region = regions.next();
+        hideStart = region[0];
+        hideEnd = region[1];
+      }
+      if (hideEnd < p)
+      {
+        hideStart = getLength();
+      }
+      // update boundary for sequence
+      if (p < hideStart)
+      {
+        start = p;
+        foundStart = true;
+      }
+    }
+
+    if (foundStart)
+    {
+      return start;
+    }
+    // otherwise, sequence was completely hidden
+    return 0;
+  }
 }