apply version 2.7 copyright
[jalview.git] / src / jalview / analysis / Dna.java
index 960a6db..49c37df 100644 (file)
@@ -1,5 +1,23 @@
+/*\r
+ * Jalview - A Sequence Alignment Editor and Viewer (Version 2.7)\r
+ * Copyright (C) 2011 J Procter, AM Waterhouse, G Barton, M Clamp, S Searle\r
+ * \r
+ * This file is part of Jalview.\r
+ * \r
+ * Jalview is free software: you can redistribute it and/or\r
+ * modify it under the terms of the GNU General Public License \r
+ * as published by the Free Software Foundation, either version 3 of the License, or (at your option) any later version.\r
+ * \r
+ * Jalview is distributed in the hope that it will be useful, but \r
+ * WITHOUT ANY WARRANTY; without even the implied warranty \r
+ * of MERCHANTABILITY or FITNESS FOR A PARTICULAR \r
+ * PURPOSE.  See the GNU General Public License for more details.\r
+ * \r
+ * You should have received a copy of the GNU General Public License along with Jalview.  If not, see <http://www.gnu.org/licenses/>.\r
+ */\r
 package jalview.analysis;\r
 \r
+import java.util.Enumeration;\r
 import java.util.Hashtable;\r
 import java.util.Vector;\r
 \r
@@ -9,6 +27,7 @@ import jalview.datamodel.AlignmentAnnotation;
 import jalview.datamodel.AlignmentI;\r
 import jalview.datamodel.Annotation;\r
 import jalview.datamodel.ColumnSelection;\r
+import jalview.datamodel.DBRefEntry;\r
 import jalview.datamodel.FeatureProperties;\r
 import jalview.datamodel.Mapping;\r
 import jalview.datamodel.Sequence;\r
@@ -34,9 +53,9 @@ public class Dna
       return 0;\r
     if (cdp1[0] < cdp2[0] || cdp1[1] < cdp2[1] || cdp1[2] < cdp2[2])\r
       return -1; // one base in cdp1 precedes the corresponding base in the\r
-                  // other codon\r
+    // other codon\r
     return 1; // one base in cdp1 appears after the corresponding base in the\r
-              // other codon.\r
+    // other codon.\r
   }\r
 \r
   /**\r
@@ -64,27 +83,60 @@ public class Dna
    * @param gapCharacter\r
    * @param annotations\r
    * @param aWidth\r
+   * @param dataset\r
+   *          destination dataset for translated sequences and mappings\r
    * @return\r
    */\r
   public static AlignmentI CdnaTranslate(SequenceI[] selection,\r
           String[] seqstring, int viscontigs[], char gapCharacter,\r
-          AlignmentAnnotation[] annotations, int aWidth)\r
+          AlignmentAnnotation[] annotations, int aWidth, Alignment dataset)\r
+  {\r
+    return CdnaTranslate(selection, seqstring, null, viscontigs,\r
+            gapCharacter, annotations, aWidth, dataset);\r
+  }\r
+\r
+  /**\r
+   * \r
+   * @param selection\r
+   * @param seqstring\r
+   * @param product\r
+   *          - array of DbRefEntry objects from which exon map in seqstring is\r
+   *          derived\r
+   * @param viscontigs\r
+   * @param gapCharacter\r
+   * @param annotations\r
+   * @param aWidth\r
+   * @param dataset\r
+   * @return\r
+   */\r
+  public static AlignmentI CdnaTranslate(SequenceI[] selection,\r
+          String[] seqstring, DBRefEntry[] product, int viscontigs[],\r
+          char gapCharacter, AlignmentAnnotation[] annotations, int aWidth,\r
+          Alignment dataset)\r
   {\r
     AlignedCodonFrame codons = new AlignedCodonFrame(aWidth); // stores hash of\r
-                                                              // subsequent\r
-                                                              // positions for\r
-                                                              // each codon\r
-                                                              // start position\r
-                                                              // in alignment\r
+    // subsequent\r
+    // positions for\r
+    // each codon\r
+    // start position\r
+    // in alignment\r
     int s, sSize = selection.length;\r
     Vector pepseqs = new Vector();\r
     for (s = 0; s < sSize; s++)\r
     {\r
       SequenceI newseq = translateCodingRegion(selection[s], seqstring[s],\r
-              viscontigs, codons, gapCharacter);\r
+              viscontigs, codons, gapCharacter,\r
+              (product != null) ? product[s] : null); // possibly anonymous\r
+      // product\r
       if (newseq != null)\r
       {\r
         pepseqs.addElement(newseq);\r
+        SequenceI ds = newseq;\r
+        while (ds.getDatasetSequence() != null)\r
+        {\r
+          ds = ds.getDatasetSequence();\r
+        }\r
+        dataset.addSequence(ds);\r
       }\r
     }\r
     if (codons.aaWidth == 0)\r
@@ -93,13 +145,139 @@ public class Dna
     pepseqs.copyInto(newseqs);\r
     AlignmentI al = new Alignment(newseqs);\r
     al.padGaps(); // ensure we look aligned.\r
-    al.setDataset(null);\r
+    al.setDataset(dataset);\r
     translateAlignedAnnotations(annotations, al, codons);\r
     al.addCodonFrame(codons);\r
     return al;\r
   }\r
 \r
   /**\r
+   * fake the collection of DbRefs with associated exon mappings to identify if\r
+   * a translation would generate distinct product in the currently selected\r
+   * region.\r
+   * \r
+   * @param selection\r
+   * @param viscontigs\r
+   * @return\r
+   */\r
+  public static boolean canTranslate(SequenceI[] selection,\r
+          int viscontigs[])\r
+  {\r
+    for (int gd = 0; gd < selection.length; gd++)\r
+    {\r
+      SequenceI dna = selection[gd];\r
+      jalview.datamodel.DBRefEntry[] dnarefs = jalview.util.DBRefUtils\r
+              .selectRefs(dna.getDBRef(),\r
+                      jalview.datamodel.DBRefSource.DNACODINGDBS);\r
+      if (dnarefs != null)\r
+      {\r
+        // intersect with pep\r
+        // intersect with pep\r
+        Vector mappedrefs = new Vector();\r
+        DBRefEntry[] refs = dna.getDBRef();\r
+        for (int d = 0; d < refs.length; d++)\r
+        {\r
+          if (refs[d].getMap() != null && refs[d].getMap().getMap() != null\r
+                  && refs[d].getMap().getMap().getFromRatio() == 3\r
+                  && refs[d].getMap().getMap().getToRatio() == 1)\r
+          {\r
+            mappedrefs.addElement(refs[d]); // add translated protein maps\r
+          }\r
+        }\r
+        dnarefs = new DBRefEntry[mappedrefs.size()];\r
+        mappedrefs.copyInto(dnarefs);\r
+        for (int d = 0; d < dnarefs.length; d++)\r
+        {\r
+          Mapping mp = dnarefs[d].getMap();\r
+          if (mp != null)\r
+          {\r
+            for (int vc = 0; vc < viscontigs.length; vc += 2)\r
+            {\r
+              int[] mpr = mp.locateMappedRange(viscontigs[vc],\r
+                      viscontigs[vc + 1]);\r
+              if (mpr != null)\r
+              {\r
+                return true;\r
+              }\r
+            }\r
+          }\r
+        }\r
+      }\r
+    }\r
+    return false;\r
+  }\r
+\r
+  /**\r
+   * generate a set of translated protein products from annotated sequenceI\r
+   * \r
+   * @param selection\r
+   * @param viscontigs\r
+   * @param gapCharacter\r
+   * @param dataset\r
+   *          destination dataset for translated sequences\r
+   * @param annotations\r
+   * @param aWidth\r
+   * @return\r
+   */\r
+  public static AlignmentI CdnaTranslate(SequenceI[] selection,\r
+          int viscontigs[], char gapCharacter, Alignment dataset)\r
+  {\r
+    int alwidth = 0;\r
+    Vector cdnasqs = new Vector();\r
+    Vector cdnasqi = new Vector();\r
+    Vector cdnaprod = new Vector();\r
+    for (int gd = 0; gd < selection.length; gd++)\r
+    {\r
+      SequenceI dna = selection[gd];\r
+      jalview.datamodel.DBRefEntry[] dnarefs = jalview.util.DBRefUtils\r
+              .selectRefs(dna.getDBRef(),\r
+                      jalview.datamodel.DBRefSource.DNACODINGDBS);\r
+      if (dnarefs != null)\r
+      {\r
+        // intersect with pep\r
+        Vector mappedrefs = new Vector();\r
+        DBRefEntry[] refs = dna.getDBRef();\r
+        for (int d = 0; d < refs.length; d++)\r
+        {\r
+          if (refs[d].getMap() != null && refs[d].getMap().getMap() != null\r
+                  && refs[d].getMap().getMap().getFromRatio() == 3\r
+                  && refs[d].getMap().getMap().getToRatio() == 1)\r
+          {\r
+            mappedrefs.addElement(refs[d]); // add translated protein maps\r
+          }\r
+        }\r
+        dnarefs = new DBRefEntry[mappedrefs.size()];\r
+        mappedrefs.copyInto(dnarefs);\r
+        for (int d = 0; d < dnarefs.length; d++)\r
+        {\r
+          Mapping mp = dnarefs[d].getMap();\r
+          StringBuffer sqstr = new StringBuffer();\r
+          if (mp != null)\r
+          {\r
+            Mapping intersect = mp.intersectVisContigs(viscontigs);\r
+            // generate seqstring for this sequence based on mapping\r
+\r
+            if (sqstr.length() > alwidth)\r
+              alwidth = sqstr.length();\r
+            cdnasqs.addElement(sqstr.toString());\r
+            cdnasqi.addElement(dna);\r
+            cdnaprod.addElement(intersect);\r
+          }\r
+        }\r
+      }\r
+      SequenceI[] cdna = new SequenceI[cdnasqs.size()];\r
+      DBRefEntry[] prods = new DBRefEntry[cdnaprod.size()];\r
+      String[] xons = new String[cdnasqs.size()];\r
+      cdnasqs.copyInto(xons);\r
+      cdnaprod.copyInto(prods);\r
+      cdnasqi.copyInto(cdna);\r
+      return CdnaTranslate(cdna, xons, prods, viscontigs, gapCharacter,\r
+              null, alwidth, dataset);\r
+    }\r
+    return null;\r
+  }\r
+\r
+  /**\r
    * translate na alignment annotations onto translated amino acid alignment al\r
    * using codon mapping codons\r
    * \r
@@ -139,18 +317,22 @@ public class Dna
             if (codons.codons[a] != null\r
                     && codons.codons[a][0] == (codons.codons[a][2] - 2))\r
             {\r
-              pos = codons.codons[a][0];\r
-              if (annotations[i].annotations[pos] == null\r
-                      || annotations[i].annotations[pos] == null)\r
-                continue;\r
-              // We just take the annotation in the first base in the codon\r
-              anots[a] = new Annotation(annotations[i].annotations[pos]);\r
+              anots[a] = getCodonAnnotation(codons.codons[a],\r
+                      annotations[i].annotations);\r
             }\r
           }\r
         }\r
 \r
         jalview.datamodel.AlignmentAnnotation aa = new jalview.datamodel.AlignmentAnnotation(\r
                 annotations[i].label, annotations[i].description, anots);\r
+        aa.graph = annotations[i].graph;\r
+        aa.graphGroup = annotations[i].graphGroup;\r
+        aa.graphHeight = annotations[i].graphHeight;\r
+        if (annotations[i].getThreshold() != null)\r
+        {\r
+          aa.setThreshold(new jalview.datamodel.GraphLine(annotations[i]\r
+                  .getThreshold()));\r
+        }\r
         if (annotations[i].hasScore)\r
         {\r
           aa.setScore(annotations[i].getScore());\r
@@ -165,7 +347,7 @@ public class Dna
             // positioning\r
             aa.setSequenceRef(aaSeq);\r
             aa.createSequenceMapping(aaSeq, aaSeq.getStart(), true); // rebuild\r
-                                                                      // mapping\r
+            // mapping\r
             aa.adjustForAlignment();\r
             aaSeq.addAlignmentAnnotation(aa);\r
           }\r
@@ -176,31 +358,92 @@ public class Dna
     }\r
   }\r
 \r
+  private static Annotation getCodonAnnotation(int[] is,\r
+          Annotation[] annotations)\r
+  {\r
+    // Have a look at all the codon positions for annotation and put the first\r
+    // one found into the translated annotation pos.\r
+    int contrib = 0;\r
+    Annotation annot = null;\r
+    for (int p = 0; p < 3; p++)\r
+    {\r
+      if (annotations[is[p]] != null)\r
+      {\r
+        if (annot == null)\r
+        {\r
+          annot = new Annotation(annotations[is[p]]);\r
+          contrib = 1;\r
+        }\r
+        else\r
+        {\r
+          // merge with last\r
+          Annotation cpy = new Annotation(annotations[is[p]]);\r
+          if (annot.colour == null)\r
+          {\r
+            annot.colour = cpy.colour;\r
+          }\r
+          if (annot.description == null || annot.description.length() == 0)\r
+          {\r
+            annot.description = cpy.description;\r
+          }\r
+          if (annot.displayCharacter == null)\r
+          {\r
+            annot.displayCharacter = cpy.displayCharacter;\r
+          }\r
+          if (annot.secondaryStructure == 0)\r
+          {\r
+            annot.secondaryStructure = cpy.secondaryStructure;\r
+          }\r
+          annot.value += cpy.value;\r
+          contrib++;\r
+        }\r
+      }\r
+    }\r
+    if (contrib > 1)\r
+    {\r
+      annot.value /= (float) contrib;\r
+    }\r
+    return annot;\r
+  }\r
+\r
   /**\r
    * Translate a na sequence\r
    * \r
    * @param selection\r
+   *          sequence displayed under viscontigs visible columns\r
    * @param seqstring\r
+   *          ORF read in some global alignment reference frame\r
    * @param viscontigs\r
+   *          mapping from global reference frame to visible seqstring ORF read\r
    * @param codons\r
+   *          Definition of global ORF alignment reference frame\r
    * @param gapCharacter\r
    * @param newSeq\r
    * @return sequence ready to be added to alignment.\r
    */\r
   public static SequenceI translateCodingRegion(SequenceI selection,\r
           String seqstring, int[] viscontigs, AlignedCodonFrame codons,\r
-          char gapCharacter)\r
+          char gapCharacter, DBRefEntry product)\r
   {\r
+    Vector skip = new Vector();\r
+    int skipint[] = null;\r
     ShiftList vismapping = new ShiftList(); // map from viscontigs to seqstring\r
-                                            // intervals\r
+    // intervals\r
     int vc, scontigs[] = new int[viscontigs.length];\r
     int npos = 0;\r
     for (vc = 0; vc < viscontigs.length; vc += 2)\r
     {\r
-      vismapping.addShift(npos, viscontigs[vc]);\r
-      scontigs[vc] = npos;\r
-      npos += viscontigs[vc + 1];\r
-      scontigs[vc + 1] = npos;\r
+      if (vc == 0)\r
+      {\r
+        vismapping.addShift(npos, viscontigs[vc]);\r
+      }\r
+      else\r
+      {\r
+        // hidden region\r
+        vismapping.addShift(npos, viscontigs[vc] - viscontigs[vc - 1] + 1);\r
+      }\r
+      scontigs[vc] = viscontigs[vc];\r
+      scontigs[vc + 1] = viscontigs[vc + 1];\r
     }\r
 \r
     StringBuffer protein = new StringBuffer();\r
@@ -222,9 +465,45 @@ public class Dna
         String aa = ResidueProperties.codonTranslate(new String(codon));\r
         rf = 0;\r
         if (aa == null)\r
+        {\r
           aa = String.valueOf(gapCharacter);\r
+          if (skipint == null)\r
+          {\r
+            skipint = new int[]\r
+            { cdp[0], cdp[2] };\r
+          }\r
+          skipint[1] = cdp[2];\r
+        }\r
         else\r
         {\r
+          if (skipint != null)\r
+          {\r
+            // edit scontigs\r
+            skipint[0] = vismapping.shift(skipint[0]);\r
+            skipint[1] = vismapping.shift(skipint[1]);\r
+            for (vc = 0; vc < scontigs.length; vc += 2)\r
+            {\r
+              if (scontigs[vc + 1] < skipint[0])\r
+              {\r
+                continue;\r
+              }\r
+              if (scontigs[vc] <= skipint[0])\r
+              {\r
+                if (skipint[0] == scontigs[vc])\r
+                {\r
+\r
+                }\r
+                else\r
+                {\r
+                  int[] t = new int[scontigs.length + 2];\r
+                  System.arraycopy(scontigs, 0, t, 0, vc - 1);\r
+                  // scontigs[vc]; //\r
+                }\r
+              }\r
+            }\r
+            skip.addElement(skipint);\r
+            skipint = null;\r
+          }\r
           if (aa.equals("STOP"))\r
           {\r
             aa = "X";\r
@@ -250,8 +529,8 @@ public class Dna
             // with a gap\r
             aa = "" + gapCharacter + aa;\r
             aspos++;\r
-            if (aspos >= codons.aaWidth)\r
-              codons.aaWidth = aspos + 1;\r
+            // if (aspos >= codons.aaWidth)\r
+            // codons.aaWidth = aspos + 1;\r
             break; // check the next position for alignment\r
           case 0:\r
             // codon aligns at aspos position.\r
@@ -267,15 +546,21 @@ public class Dna
           codons.codons[aspos] = new int[]\r
           { cdp[0], cdp[1], cdp[2] };\r
         }\r
-        aspos++;\r
         if (aspos >= codons.aaWidth)\r
-          codons.aaWidth = aspos + 1;\r
+        {\r
+          // update maximum alignment width\r
+          // (we can do this without calling checkCodonFrameWidth because it was\r
+          // already done above)\r
+          codons.setAaWidth(aspos);\r
+        }\r
+        // ready for next translated reading frame alignment position (if any)\r
+        aspos++;\r
       }\r
     }\r
     if (resSize > 0)\r
     {\r
-      SequenceI newseq = new Sequence(selection.getName(), protein\r
-              .toString());\r
+      SequenceI newseq = new Sequence(selection.getName(),\r
+              protein.toString());\r
       if (rf != 0)\r
       {\r
         jalview.bin.Cache.log\r
@@ -283,8 +568,8 @@ public class Dna
         // map and trim contigs to ORF region\r
         vc = scontigs.length - 1;\r
         lastnpos = vismapping.shift(lastnpos); // place npos in context of\r
-                                                // whole dna alignment (rather\r
-                                                // than visible contigs)\r
+        // whole dna alignment (rather\r
+        // than visible contigs)\r
         // incomplete ORF could be broken over one or two visible contig\r
         // intervals.\r
         while (vc >= 0 && scontigs[vc] > lastnpos)\r
@@ -313,14 +598,11 @@ public class Dna
       if (scontigs != null)\r
       {\r
         npos = 0;\r
-        // Find sequence position for scontigs positions on the nucleotide\r
-        // sequence string we were passed.\r
-        for (vc = 0; vc < viscontigs.length; vc += 2)\r
+        // map scontigs to actual sequence positions on selection\r
+        for (vc = 0; vc < scontigs.length; vc += 2)\r
         {\r
           scontigs[vc] = selection.findPosition(scontigs[vc]); // not from 1!\r
-          npos += viscontigs[vc];\r
-          scontigs[vc + 1] = selection\r
-                  .findPosition(npos + scontigs[vc + 1]); // exclusive\r
+          scontigs[vc + 1] = selection.findPosition(scontigs[vc + 1]); // exclusive\r
           if (scontigs[vc + 1] == selection.getEnd())\r
             break;\r
         }\r
@@ -331,22 +613,49 @@ public class Dna
           System.arraycopy(scontigs, 0, t, 0, vc + 2);\r
           scontigs = t;\r
         }\r
-\r
+        /*\r
+         * delete intervals in scontigs which are not translated. 1. map skip\r
+         * into sequence position intervals 2. truncate existing ranges and add\r
+         * new ranges to exclude untranslated regions. if (skip.size()>0) {\r
+         * Vector narange = new Vector(); for (vc=0; vc<scontigs.length; vc++) {\r
+         * narange.addElement(new int[] {scontigs[vc]}); } int sint=0,iv[]; vc =\r
+         * 0; while (sint<skip.size()) { skipint = (int[]) skip.elementAt(sint);\r
+         * do { iv = (int[]) narange.elementAt(vc); if (iv[0]>=skipint[0] &&\r
+         * iv[0]<=skipint[1]) { if (iv[0]==skipint[0]) { // delete beginning of\r
+         * range } else { // truncate range and create new one if necessary iv =\r
+         * (int[]) narange.elementAt(vc+1); if (iv[0]<=skipint[1]) { // truncate\r
+         * range iv[0] = skipint[1]; } else { } } } else if (iv[0]<skipint[0]) {\r
+         * iv = (int[]) narange.elementAt(vc+1); } } while (iv[0]) } }\r
+         */\r
         MapList map = new MapList(scontigs, new int[]\r
-        { 1, resSize }, 3, 1); // TODO: store mapping on newSeq for linked\r
-                                // DNA/Protein viewing.\r
+        { 1, resSize }, 3, 1);\r
+\r
+        // update newseq as if it was generated as mapping from product\r
+\r
+        if (product != null)\r
+        {\r
+          newseq.setName(product.getSource() + "|"\r
+                  + product.getAccessionId());\r
+          if (product.getMap() != null)\r
+          {\r
+            // Mapping mp = product.getMap();\r
+            // newseq.setStart(mp.getPosition(scontigs[0]));\r
+            // newseq.setEnd(mp\r
+            // .getPosition(scontigs[scontigs.length - 1]));\r
+          }\r
+        }\r
         transferCodedFeatures(selection, newseq, map, null, null);\r
         SequenceI rseq = newseq.deriveSequence(); // construct a dataset\r
-                                                  // sequence for our new\r
-                                                  // peptide, regardless.\r
+        // sequence for our new\r
+        // peptide, regardless.\r
         // store a mapping (this actually stores a mapping between the dataset\r
         // sequences for the two sequences\r
-        codons.addMap(selection, newseq, map);\r
+        codons.addMap(selection, rseq, map);\r
         return rseq;\r
       }\r
     }\r
     // register the mapping somehow\r
-    // \r
+    //\r
     return null;\r
   }\r
 \r
@@ -391,8 +700,7 @@ public class Dna
         fgstate = (featureGroups == null) ? null : ((Boolean) featureGroups\r
                 .get(sf[f].featureGroup));\r
         if ((featureTypes == null || featureTypes.containsKey(sf[f]\r
-                .getType()))\r
-                && (fgstate == null || fgstate.booleanValue()))\r
+                .getType())) && (fgstate == null || fgstate.booleanValue()))\r
         {\r
           if (FeatureProperties.isCodingFeature(null, sf[f].getType()))\r
           {\r